SET UP

# load packages
library(tidyverse) #data wrangling
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.0     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.1     ✔ tibble    3.2.0
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.0
## ✔ purrr     1.0.1     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the ]8;;http://conflicted.r-lib.org/conflicted package]8;; to force all conflicts to become errors
library(dplyr) # data wrangling
library(ggpubr) # plotting
library(readr) # read_csv
library(psych) # psych statistics - MSSD
## 
## Attaching package: 'psych'
## 
## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
library(lme4) # mixed models
## Loading required package: Matrix
## 
## Attaching package: 'Matrix'
## 
## The following objects are masked from 'package:tidyr':
## 
##     expand, pack, unpack
library(effectsize) #effect size functions
## 
## Attaching package: 'effectsize'
## 
## The following object is masked from 'package:psych':
## 
##     phi
library(EMAtools) #effect size functions
library(lubridate) #time variables
library(tseries) #time series tools
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

EMA

filename <- "20230321_EMAData.csv"

EMAdata <- read_csv(filename)
## Rows: 1562 Columns: 34
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr   (1): date
## dbl  (30): id, event, day, instance, enrollmentgroup, completed, positive_sc...
## time  (3): time, time_participant, time_interval
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
filename <- "20230327_betweenday.csv"

betweenday <- read_csv(filename)
## Rows: 286 Columns: 16
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr   (1): date
## dbl  (12): id, event, day, instance, enrollmentgroup, completed, positive_sc...
## time  (3): time, time_participant, time_interval
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
#transform variables to factor
EMAdata$id <- as.factor(EMAdata$id)
EMAdata$event <- as.factor(EMAdata$event) 
levels(EMAdata$event) <- c("Baseline", "Session 4", "Post")
EMAdata$enrollmentgroup <- as.factor(EMAdata$enrollmentgroup) 
levels(EMAdata$enrollmentgroup) <- c("Waitlist Control", "Immediate Enrollment", "Group")

Descriptive Statistics

#remove group 2
EMAdata <- filter(EMAdata, enrollmentgroup != "Group")

Missing Data

Missing Instances

#create dataframe of missing instances (completed = 0)
missingdata <- filter(EMAdata, completed == '0')
#create table of n missing: event x group
missing_instance <- matrix(0, 4,4)

colnames(missing_instance) <- c("Baseline", "Session4", "Post", "Total")
rownames(missing_instance) <- c("WLC", "IE", "Group", "Total")
missing_instance[1,1] <- length(which(missingdata$event == "baseline" & missingdata$enrollmentgroup == "0"))
missing_instance[2,1] <- length(which(missingdata$event == "baseline" & missingdata$enrollmentgroup == "1"))
missing_instance[3,1] <- length(which(missingdata$event == "baseline" & missingdata$enrollmentgroup == "2"))
missing_instance[1,2] <- length(which(missingdata$event == "session4" & missingdata$enrollmentgroup == "0"))
missing_instance[2,2] <- length(which(missingdata$event == "session4" & missingdata$enrollmentgroup == "1"))
missing_instance[3,2] <- length(which(missingdata$event == "session4" & missingdata$enrollmentgroup == "2"))
missing_instance[1,3] <- length(which(missingdata$event == "post" & missingdata$enrollmentgroup == "0"))
missing_instance[2,3] <- length(which(missingdata$event == "post" & missingdata$enrollmentgroup == "1"))
missing_instance[3,3] <- length(which(missingdata$event == "post" & missingdata$enrollmentgroup == "2"))

missing_instance[1,4] <- length(which(missingdata$enrollmentgroup == "0"))
missing_instance[2,4] <- length(which(missingdata$enrollmentgroup == "1"))
missing_instance[3,4] <- length(which(missingdata$enrollmentgroup == "2"))

missing_instance[4,1] <- length(which(missingdata$event == "baseline"))
missing_instance[4,2] <- length(which(missingdata$event == "session4"))
missing_instance[4,3] <- length(which(missingdata$event == "post"))
missing_instance[4,4] <- missing_instance[4,1] + missing_instance[4,2] + missing_instance[4,3]

Missing Items

#create dataframe of missing items (completed =1, N/A on at least 1 item)
missingitems <- filter(EMAdata, completed == "1" & (is.na(positive_score) | is.na(negative_score))) #n=9
#impute missing items using the average of other neg/pos items
## create function

#function inputs: ID, EVENT, INSTANCE
itemavg <- function(ID, EVENT, INSTANCE){

  # row from id, event, instance
  r <- which(EMAdata$id == ID & EMAdata$event == EVENT & EMAdata$instance == INSTANCE)  
  
  positems <- c()
  negitems <- c()

  #positive items: 1, 3, 5, 9, 10, 12, 14, 16, 17, 19
  positems <- append(positems, as.integer(EMAdata[r,15])) #posq1
  positems <- append(positems, as.integer(EMAdata[r,17])) #posq3
  positems <- append(positems, as.integer(EMAdata[r,19]))#posq5
  positems <- append(positems, as.integer(EMAdata[r,23]))#posq9
  positems <- append(positems, as.integer(EMAdata[r,24])) #posq10
  positems <- append(positems, as.integer(EMAdata[r,26])) #posq12
  positems <- append(positems, as.integer(EMAdata[r,28])) #posq14
  positems <- append(positems, as.integer(EMAdata[r,30])) #posq16
  positems <- append(positems, as.integer(EMAdata[r,31])) #posq17
  positems <- append(positems, as.integer(EMAdata[r,33])) #posq19

  posavg <- mean(positems, na.rm=TRUE)

  #negative items: 2, 4, 6, 7, 8, 11, 13, 15, 18, 20
  negitems <- append(negitems,  as.integer(EMAdata[r,16])) #negq2
  negitems <- append(negitems, as.integer(EMAdata[r,18])) #negq4
  negitems <- append(negitems, as.integer(EMAdata[r,20])) #negq6
  negitems <- append(negitems, as.integer(EMAdata[r,21])) #negq7
  negitems <- append(negitems, as.integer(EMAdata[r,22])) #negq8
  negitems <- append(negitems, as.integer(EMAdata[r,25])) #negq11
  negitems <- append(negitems, as.integer(EMAdata[r,27])) #negq13
  negitems <- append(negitems, as.integer(EMAdata[r,29])) #negq15
  negitems <- append(negitems, as.integer(EMAdata[r,32])) #negq18
  negitems <- append(negitems, as.integer(EMAdata[r,34])) #negq20
  
  negavg <- mean(negitems, na.rm=TRUE)
  
  #output posavg, negavg
  out <- list(r, posavg, negavg)
  
  return(out)
}  
#add all positive items
posscore <- function(r){
  score <- as.integer(EMAdata[r, 15]) + as.integer(EMAdata[r, 17]) + 
    as.integer(EMAdata[r, 19]) + as.integer(EMAdata[r, 23]) +
    as.integer(EMAdata[r, 24]) + as.integer(EMAdata[r, 26]) + 
    as.integer(EMAdata[r, 28]) + as.integer(EMAdata[r, 30]) + 
    as.integer(EMAdata[r, 31]) + as.integer(EMAdata[r, 33])
  
  return(score)
}

#add all negative items  
negscore <- function(r){
  score <- as.integer(EMAdata[r, 16]) + as.integer(EMAdata[r, 18]) + 
    as.integer(EMAdata[r, 20]) + as.integer(EMAdata[r, 21]) +
    as.integer(EMAdata[r, 22]) + as.integer(EMAdata[r, 25]) + 
    as.integer(EMAdata[r, 27]) + as.integer(EMAdata[r, 29]) + 
    as.integer(EMAdata[r, 32]) + as.integer(EMAdata[r, 34])
  
  return(score)
}
#12015 session 4:8
out<- itemavg("12015", "Session 4", 8)
EMAdata[23, 17] <- out[2] #posq3 
EMAdata[23, 23] <- out[2] #posq9 
EMAdata[23,11] <- posscore(23) #positive score

#12015 session4: 12
out<- itemavg("12015", "Session 4", 12)
EMAdata[27, 25] <- out[3] #negq11
EMAdata[27,12] <- negscore(27) #negative score

#12017 post: 3
out <- itemavg("12017", "Post", 3)
EMAdata[123, 30] <- out[2] #posq16
EMAdata[123, 11] <- posscore(123)

#12027 baseline:2
out <- itemavg("12027", "Baseline", 2)
EMAdata[542, 28] <- out[2] #posq14
EMAdata[542, 11] <- posscore(542)

#12027 session4:12
out <- itemavg("12027", "Session 4", 12)
EMAdata[567, 22] <- out[3] #neg8
EMAdata[567, 27] <- out[3] #neg13
EMAdata[567,12] <- negscore(567)

#12033 baseline:15
out<- itemavg("12033", "Baseline", 15)
EMAdata[810, 24] <- out[2] #pos10
EMAdata[810, 11] <- posscore(810)

#12033 session4:11
out<- itemavg("12033", "Session 4", 11)
EMAdata[821, 30] <- out[2] #posq16
EMAdata[821, 11] <- posscore(821)

#12039 post:2
out<- itemavg("12039", "Post", 2)
EMAdata[1097, 31] <- out[2] #posq17
EMAdata[1097,11] <- posscore(1097)

#12039 post:3
out<- itemavg("12039", "Post", 3)
EMAdata[1098, 28] <- out[2] #posq14
EMAdata[1098, 34] <- out[3] #neg20
EMAdata[1098, 11] <- posscore(1098)
EMAdata[1098, 12] <- negscore(1098)
#check missing items (completed =1, N/A on at least 1 item)
missingitems <- filter(EMAdata, completed == "1" & (is.na(positive_score) | is.na(negative_score))) #n=9
#save to csv
write.csv(EMAdata, "20230411_EMAData.csv")

Survey Completion

Compliance Rate

#create dataframe of id x completion
compliance_rate <- matrix(0, 35, 3)

colnames(compliance_rate) <- c("Group", "Completed", "ComplianceRate")
id_names <- c("12015", "12016", "12017", "12018", "12019", "12020","12021", "12022", "12023", "12024","12025", "12026", "12027", "12028", "12029", "12030", "12031", "12032", "12033", "12034", "12035", "12036", "12037", "12038", "12039", "12040", "12041","12043", "12044", "12045", "12048","12050", "12055","12056", "12057")
rownames(compliance_rate) <- id_names
#fill matrix col 1 with enrollment group
r = 1
n = 1

for(n in 1:35){
  ind <- which(EMAdata$id == id_names[n])
  compliance_rate[r,1] <- EMAdata$enrollmentgroup[ind[1]]
  r <- r + 1
  n <- n + 1
}
#fill matrix col 2 with n completed by id
r = 1
n = 1

for(n in 1:35){
  compliance_rate[r,2] <- length(which(EMAdata$id == id_names[n]))
  r <- r + 1
  n <- n + 1
  }
#fill matrix col 3 with compliance rate
r = 1
n = 1

for (n in 1:35){
  compliance_rate[r, 3] <- compliance_rate[r, 2]/45
  r <- r + 1
  n <- n + 1
}
#output min and max values
max <- which(compliance_rate == max(compliance_rate[,3]), arr.ind=TRUE)
min <- which(compliance_rate == min(compliance_rate[,3]), arr.ind=TRUE)

print(paste("Maximum compliance rate:", compliance_rate[max[1]]))
## [1] "Maximum compliance rate: 1"
print(paste("Minimum compliance rate:", compliance_rate[min]))
## [1] "Minimum compliance rate: 0" "Minimum compliance rate: 0"
## [3] "Minimum compliance rate: 0" "Minimum compliance rate: 0"
## [5] "Minimum compliance rate: 0" "Minimum compliance rate: 0"

Mean (SD) Table

summarise_at(EMAdata, vars(positive_score, negative_score), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 1 × 4
##   positive_score_mean negative_score_mean positive_score_sd negative_score_sd
##                 <dbl>               <dbl>             <dbl>             <dbl>
## 1                22.5                15.9              7.82              4.61
#group by event
EMAdata%>%
  group_by(event) %>%
  summarise_at(vars(positive_score, negative_score), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 3 × 5
##   event     positive_score_mean negative_score_mean positive_score_sd negative…¹
##   <fct>                   <dbl>               <dbl>             <dbl>      <dbl>
## 1 Baseline                 21.5                16.0              7.10       4.73
## 2 Session 4                21.8                16.0              7.56       4.89
## 3 Post                     24.4                15.7              8.43       4.15
## # … with abbreviated variable name ¹​negative_score_sd
#group by enrollment
EMAdata%>%
  group_by(enrollmentgroup) %>%
  summarise_at(vars(positive_score, negative_score), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 2 × 5
##   enrollmentgroup      positive_score_mean negative_score_mean positiv…¹ negat…²
##   <fct>                              <dbl>               <dbl>     <dbl>   <dbl>
## 1 Waitlist Control                    22.8                16.0      7.97    4.23
## 2 Immediate Enrollment                22.3                15.8      7.67    4.94
## # … with abbreviated variable names ¹​positive_score_sd, ²​negative_score_sd

Normality & Distribution

#create grouped boxplot - positive

ggplot(EMAdata, aes(x=event, y=positive_score, fill=enrollmentgroup)) + scale_fill_manual(values = c("lightblue3", "dodgerblue4")) + geom_boxplot() + xlab("event") + ylab("score") + ggtitle("PANAS Positive Score x Group")
## Warning: Removed 185 rows containing non-finite values (`stat_boxplot()`).

#create grouped boxplot - negative

ggplot(EMAdata, aes(x=event, y=negative_score, fill=enrollmentgroup)) + 
  scale_fill_manual(values = c("lightblue3", "dodgerblue4")) +
  geom_boxplot() + xlab("event") + ylab("score") +
  ggtitle("PANAS Negative Score x Group")
## Warning: Removed 185 rows containing non-finite values (`stat_boxplot()`).

## Affect Dynamic Variables

withinday <- EMAdata[,c('id', 'event', 'day', 'instance', 'enrollmentgroup', 'time_interval', 'positive_score', 'negative_score', "positive_SD","negative_SD")]
#convert time_interval from h:m:s to minutes
withinday$time_interval_min <- hour(withinday$time_interval)*60 + minute(withinday$time_interval)

withinday <- withinday %>%
  relocate("time_interval_min", .after="time_interval")
# scatter plot of ASD x time interval (hour)
ggplot(withinday, aes(x = time_interval, y=positive_SD)) + geom_point(size = 0.5) +
  ggtitle("Successive Difference - Positive")
## Warning: Removed 406 rows containing missing values (`geom_point()`).

ggplot(withinday, aes(x = time_interval, y=negative_SD)) + geom_point(size = 0.5) +
  ggtitle("Successive Difference - Negative")
## Warning: Removed 406 rows containing missing values (`geom_point()`).

#ASD = positive_SD[i] / time_interval[i]
pos_ASD <- vector()

for(i in 1:1427){
  pos_ASD <- append(pos_ASD, c(withinday$positive_SD[i] / withinday$time_interval_min[i]))
  
  i = i + 1
}

head(pos_ASD)
## [1]          NA 0.007246377 0.000000000 0.003745318 0.024752475          NA
#create positive ASD variable
withinday$positive_ASD <- pos_ASD
neg_ASD <- vector()

for(i in 1:1427){
  neg_ASD <- append(neg_ASD, c(withinday$negative_SD[i] / withinday$time_interval_min[i]))
  
  i = i + 1
}

head(neg_ASD)
## [1]         NA 0.01449275 0.00000000 0.00000000 0.00000000         NA
#create negative ASD variable
withinday$negative_ASD <- neg_ASD
# scatter plot of ASD x time interval (hour)
ggplot(withinday, aes(x = time_interval, y=positive_ASD)) + geom_point(size = 0.5) +
  ggtitle("Adj. Successive Difference - Positive")
## Warning: Removed 413 rows containing missing values (`geom_point()`).

ggplot(withinday, aes(x = time_interval, y=negative_ASD)) + geom_point(size = 0.5) +
  ggtitle("Adj. Successive Difference - Negative")
## Warning: Removed 413 rows containing missing values (`geom_point()`).

betweenday <- filter(betweenday, betweenday$id != "12046")
  
#split betweenday df by event
betweenday_baseline <- filter(betweenday, betweenday$event == "1")
betweenday_session4 <- filter(betweenday, betweenday$event == "2")
betweenday_post <- filter(betweenday, betweenday$event == "3")
#remove group 2 & 12046 (no between data)
id_names <- c("12015", "12016", "12017", "12018", "12019", "12020","12021", "12022", "12023", "12024","12025", "12026", "12027", "12028", "12029", "12030", "12031", "12032", "12033", "12034", "12035", "12036", "12037", "12038", "12039", "12040", "12041","12043", "12044", "12045", "12050", "12057")

length(id_names)
## [1] 32
id_group <- c("1", "1", "0", "1", "0", "1", "0", "0", "1", "1", "0", "1", "0", "0", "0", "0", "1", "1", "1", "0", "1", "0", "0", "1", "1", "0", "0", "1", "0", "1","1", "1")

length(id_group)
## [1] 32

Mean

# initialize MEAN matrix
ema_mean <- matrix(0, 32, 7)

colnames(ema_mean) <- c("group", "posbaseline", "possession4", "pospost", "negbaseline", "negsession4", "negpost")
rownames(ema_mean) <- id_names

for(i in 1:32){
  ema_mean[i,1] <- id_group[i]
}
#fill positive mean
## baseline

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_baseline$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_baseline$positive_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_mean[i,2] <- mean(vec)
}

## session4

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_session4$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_session4$positive_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_mean[i,3] <- mean(vec)
}

## post

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_post$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_post$positive_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_mean[i,4] <- mean(vec)
}
#fill negative mean
## baseline

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_baseline$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_baseline$negative_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_mean[i,5] <- mean(vec)
}

## session4

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_session4$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_session4$negative_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_mean[i,6] <- mean(vec)
}

## post

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_post$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_post$negative_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_mean[i,7] <- mean(vec)
}
ema_mean <- as.data.frame(ema_mean)
ema_mean$posbaseline <- as.numeric(ema_mean$posbaseline)
ema_mean$possession4 <- as.numeric(ema_mean$possession4)
ema_mean$pospost <- as.numeric(ema_mean$pospost)
ema_mean$negbaseline <- as.numeric(ema_mean$negbaseline)
ema_mean$negsession4 <- as.numeric(ema_mean$negsession4)
ema_mean$negpost <- as.numeric(ema_mean$negpost)
#group by group
ema_mean%>%
  group_by(group) %>%
  summarise_at(vars(posbaseline, possession4, pospost, negbaseline, negsession4, negpost, ), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 2 × 13
##   group posbas…¹ posse…² pospo…³ negba…⁴ negse…⁵ negpo…⁶ posba…⁷ posse…⁸ pospo…⁹
##   <chr>    <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>
## 1 0         21.7    23.3    23.2    16.5    15.9    16.6    5.61    6.10    6.16
## 2 1         22.3    20.7    23.9    17.3    16.0    16.2    5.77    5.30    5.87
## # … with 3 more variables: negbaseline_sd <dbl>, negsession4_sd <dbl>,
## #   negpost_sd <dbl>, and abbreviated variable names ¹​posbaseline_mean,
## #   ²​possession4_mean, ³​pospost_mean, ⁴​negbaseline_mean, ⁵​negsession4_mean,
## #   ⁶​negpost_mean, ⁷​posbaseline_sd, ⁸​possession4_sd, ⁹​pospost_sd

Plots

#create id variable
ema_mean$id <- id_names
  
#reformat data: baseline
ema_mean_baseline <- ema_mean %>%
  dplyr::select(id, group, pos = posbaseline, neg = negbaseline)

#create time variable
ema_mean_baseline$time <- as.integer(1)
ema_mean_baseline <- relocate(ema_mean_baseline, time, .before=group)
#reformat data: session4
ema_mean_session4 <- ema_mean %>%
  dplyr::select(id, group, pos = possession4, neg = negsession4)

#create time variable
ema_mean_session4$time <- as.integer(2)
ema_mean_session4 <- relocate(ema_mean_session4, time, .before=group)
#reformat data: post
ema_mean_post <- ema_mean %>%
  dplyr::select(id, group, pos = pospost, neg = negpost)

#create time variable
ema_mean_post$time <- as.integer(3)
ema_mean_post<- relocate(ema_mean_post, time, .before=group)
#join predata and postdata
ema_mean2 <- merge(ema_mean_baseline, ema_mean_session4, all=TRUE)
ema_mean2 <- merge(ema_mean2, ema_mean_post, all=TRUE)

ema_mean2$id <- factor(ema_mean2$id)
ema_mean2$group <- factor(ema_mean2$group)
levels(ema_mean2$group) <- c("Control", "Treatment")
ema_mean2$time <- factor(ema_mean2$time)
levels(ema_mean2$time) <- c("baseline", "session4", "post")
#plot positive means
meanpos_ggplot<- ggplot(data= ema_mean2, 
       aes(y=pos, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Mean Positive Affect Across Timepoints") + 
  labs(y = "positive affect score", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
meanpos_ggplot
## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Warning: Removed 1 row containing missing values (`geom_line()`).

#plot negative means
meanneg_ggplot<- ggplot(data= ema_mean2, 
       aes(y=neg, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Mean Negative Affect Across Timepoints") + 
  labs(y = "negative affect score", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
meanneg_ggplot
## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Warning: Removed 1 row containing missing values (`geom_line()`).

Variability (SD)

# initialize SD matrix
ema_sd <- matrix(0, 32, 7)

colnames(ema_sd) <- c("group", "posbaseline", "possession4", "pospost", "negbaseline", "negsession4", "negpost")
rownames(ema_sd) <- id_names

for(i in 1:32){
  ema_sd[i,1] <- id_group[i]
}
#fill positive mean
## baseline

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_baseline$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_baseline$positive_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_sd[i,2] <- sd(vec)
}

## session4

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_session4$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_session4$positive_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_sd[i,3] <- sd(vec)
}

## post

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_post$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_post$positive_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_sd[i,4] <- sd(vec)
}
#fill negative mean
## baseline

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_baseline$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_baseline$negative_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_sd[i,5] <- sd(vec)
}

## session4

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_session4$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_session4$negative_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_sd[i,6] <- sd(vec)
}

## post

for(i in 1:32){
  name <- id_names[i]
  index <- which(betweenday_post$id == name)
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){
    n <- betweenday_post$negative_mean[index[x]]
    vec <- c(vec, n)
  }
  
  ema_sd[i,7] <- sd(vec)
}
ema_sd <- as.data.frame(ema_sd)
ema_sd$posbaseline <- as.numeric(ema_sd$posbaseline)
ema_sd$possession4 <- as.numeric(ema_sd$possession4)
ema_sd$pospost <- as.numeric(ema_sd$pospost)
ema_sd$negbaseline <- as.numeric(ema_sd$negbaseline)
ema_sd$negsession4 <- as.numeric(ema_sd$negsession4)
ema_sd$negpost <- as.numeric(ema_sd$negpost)
#group by group
ema_sd%>%
  group_by(group) %>%
  summarise_at(vars(posbaseline, possession4, pospost, negbaseline, negsession4, negpost, ), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 2 × 13
##   group posbas…¹ posse…² pospo…³ negba…⁴ negse…⁵ negpo…⁶ posba…⁷ posse…⁸ pospo…⁹
##   <chr>    <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>
## 1 0         3.25    3.79    2.87    2.40    2.00    2.45    1.87    2.31    2.10
## 2 1         3.47    3.94    3.01    2.59    1.95    1.84    2.42    2.60    1.44
## # … with 3 more variables: negbaseline_sd <dbl>, negsession4_sd <dbl>,
## #   negpost_sd <dbl>, and abbreviated variable names ¹​posbaseline_mean,
## #   ²​possession4_mean, ³​pospost_mean, ⁴​negbaseline_mean, ⁵​negsession4_mean,
## #   ⁶​negpost_mean, ⁷​posbaseline_sd, ⁸​possession4_sd, ⁹​pospost_sd

Plots

#create id variable
ema_sd$id <- id_names
  
#reformat data: baseline
ema_sd_baseline <- ema_sd %>%
  dplyr::select(id, group, pos = posbaseline, neg = negbaseline)

#create time variable
ema_sd_baseline$time <- as.integer(1)
ema_sd_baseline <- relocate(ema_sd_baseline, time, .before=group)
#reformat data: session4
ema_sd_session4 <- ema_mean %>%
  dplyr::select(id, group, pos = possession4, neg = negsession4)

#create time variable
ema_sd_session4$time <- as.integer(2)
ema_sd_session4 <- relocate(ema_sd_session4, time, .before=group)
#reformat data: post
ema_sd_post <- ema_mean %>%
  dplyr::select(id, group, pos = pospost, neg = negpost)

#create time variable
ema_sd_post$time <- as.integer(3)
ema_sd_post<- relocate(ema_sd_post, time, .before=group)
#join predata and postdata
ema_sd2 <- merge(ema_sd_baseline, ema_sd_session4, all=TRUE)
ema_sd2  <- merge(ema_sd2 , ema_sd_post, all=TRUE)

ema_sd2 $id <- factor(ema_sd2 $id)
ema_sd2 $group <- factor(ema_sd2 $group)
levels(ema_sd2 $group) <- c("Control", "Treatment")
ema_sd2 $time <- factor(ema_sd2 $time)
levels(ema_sd2 $time) <- c("baseline", "session4", "post")
#plot positive means
sdpos_ggplot<- ggplot(data= ema_sd2, 
       aes(y=pos, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("SD of Positive Affect Across Timepoints") + 
  labs(y = "SD", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
sdpos_ggplot
## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Warning: Removed 1 row containing missing values (`geom_line()`).

#plot negative means
sdneg_ggplot<- ggplot(data= ema_sd2, 
       aes(y=neg, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("SD of Negative Affect Across Timepoints") + 
  labs(y = "SD", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
sdneg_ggplot
## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Warning: Removed 1 row containing missing values (`geom_line()`).

### Intertia (Auto-Correlation)

# create AC df
EMAdata_ac <- EMAdata[,1:14]
EMAdata_ac$positive_ac <- c(0)
EMAdata_ac$negative_ac <- c(0)

EMAdata_ac <- filter(EMAdata_ac, !is.na(positive_score))
#compute AC for positive score
for(n in 1:32){ #cycle through participant ids
  name <- id_names[n]
  
  for(d in 1:9){ #cycle through days
    index <- which(EMAdata_ac$id == name & EMAdata_ac$day == d)
      
    vec <- c()
    len <- length(index)
  
    if(len !=0){#issue if there is no entry for a day number
      for(i in 1:len){ #number of instances
        n <- EMAdata_ac$positive_score[index[i]]
        vec <- c(vec,n) #create vector of pos scores for each instance
      
      }
      acf <- unlist(acf(vec, pl=FALSE))
    
      for(i in 1:len){ #number of instances
        EMAdata_ac$positive_ac[index[i]] <- acf[i]
      }
    }
  }
}
#compute AC for negative score
for(n in 1:32){ #cycle through participant ids
  name <- id_names[n]
  
  for(d in 1:9){ #cycle through days
    index <- which(EMAdata_ac$id == name & EMAdata_ac$day == d)
      
    vec <- c()
    len <- length(index)
  
    if(len !=0){#issue if there is no entry for a day number
      for(i in 1:len){ #number of instances
        n <- EMAdata_ac$negative_score[index[i]]
        vec <- c(vec,n) #create vector of neg scores for each instance
      
      }
      acf <- unlist(acf(vec, pl=FALSE))
    
      for(i in 1:len){ #number of instances
        EMAdata_ac$negative_ac[index[i]] <- acf[i]
      }
    }
  }
}
# initialize AC matrix
ema_ac <- matrix(0, 32, 7)

colnames(ema_ac) <- c("group", "posbaseline", "possession4", "pospost", "negbaseline", "negsession4", "negpost")
rownames(ema_ac) <- id_names

for(i in 1:32){
  ema_ac[i,1] <- id_group[i]
}
EMAdata_ac$positive_ac <- as.numeric(EMAdata_ac$positive_ac)
EMAdata_ac$negative_ac <- as.numeric(EMAdata_ac$negative_ac)
#fill positive AC
## baseline

for(i in 1:32){
  name <- id_names[i]
  index <- which(EMAdata_ac$id == name & EMAdata_ac$event == "Baseline")
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){ 
    n <- EMAdata_ac$positive_ac[index[x]]
   
    if(n != "1"){
      vec <- c(vec, n)}
  }

  ema_ac[i,2] <- mean(vec, na.rm=TRUE)
}
## session4
for(i in 1:32){
  name <- id_names[i]
  index <- which(EMAdata_ac$id == name & EMAdata_ac$event == "Session 4")
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){ 
    n <- EMAdata_ac$positive_ac[index[x]]
   
    if(n != "1"){
      vec <- c(vec, n)
    }
  }
ema_ac[i,3] <- mean(vec)
}
## post
for(i in 1:32){
  name <- id_names[i]
  index <- which(EMAdata_ac$id == name & EMAdata_ac$event == "Post")
  
  vec <- c()
  len <- length(index)
  
  if(len !=0){
    for(x in 1:len){ 
      n <- EMAdata_ac$positive_ac[index[x]]
   
      if(n != "1"){
        vec <- c(vec, n)
      }
    }
  }
ema_ac[i,4] <- mean(vec, na.rm=TRUE)
}
## Warning in mean.default(vec, na.rm = TRUE): argument is not numeric or logical:
## returning NA
#fill negative AC
## baseline

for(i in 1:32){
  name <- id_names[i]
  index <- which(EMAdata_ac$id == name & EMAdata_ac$event == "Baseline")
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){ 
    n <- EMAdata_ac$negative_ac[index[x]]
   
    if(n != "1"){
      vec <- c(vec, n)}
  }

  ema_ac[i,5] <- mean(vec, na.rm=TRUE)
}
## session4

for(i in 1:32){
  name <- id_names[i]
  index <- which(EMAdata_ac$id == name & EMAdata_ac$event == "Session 4")
  
  vec <- c()
  len <- length(index)

  for(x in 1:len){ 
    n <- EMAdata_ac$negative_ac[index[x]]
   
    if(n != "1"){
      vec <- c(vec, n)
    }
  }
ema_ac[i,6] <- mean(vec)
}
## post
for(i in 1:32){
  name <- id_names[i]
  index <- which(EMAdata_ac$id == name & EMAdata_ac$event == "Post")
  
  vec <- c()
  len <- length(index)
  
  if(len !=0){
    for(x in 1:len){ 
      n <- EMAdata_ac$negative_ac[index[x]]
   
      if(n != "1"){
        vec <- c(vec, n)
      }
    }
  }
ema_ac[i,7] <- mean(vec, na.rm=TRUE)
}
## Warning in mean.default(vec, na.rm = TRUE): argument is not numeric or logical:
## returning NA
ema_ac <- as.data.frame(ema_ac)
ema_ac$posbaseline <- as.numeric(ema_ac$posbaseline)
ema_ac$possession4 <- as.numeric(ema_ac$possession4)
ema_ac$pospost <- as.numeric(ema_ac$pospost)
ema_ac$negbaseline <- as.numeric(ema_ac$negbaseline)
ema_ac$negsession4 <- as.numeric(ema_ac$negsession4)
ema_ac$negpost <- as.numeric(ema_ac$negpost)
#group by group
ema_ac%>%
  group_by(group) %>%
  summarise_at(vars(posbaseline, possession4, pospost, negbaseline, negsession4, negpost, ), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 2 × 13
##   group posbas…¹ posse…² pospo…³ negba…⁴ negse…⁵ negpo…⁶ posba…⁷ posse…⁸ pospo…⁹
##   <chr>    <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>
## 1 0       -0.152  -0.140  -0.143  -0.152  -0.140  -0.143  0.0345  0.0218  0.0258
## 2 1       -0.172  -0.149  -0.189  -0.174  -0.149  -0.189  0.0882  0.0560  0.124 
## # … with 3 more variables: negbaseline_sd <dbl>, negsession4_sd <dbl>,
## #   negpost_sd <dbl>, and abbreviated variable names ¹​posbaseline_mean,
## #   ²​possession4_mean, ³​pospost_mean, ⁴​negbaseline_mean, ⁵​negsession4_mean,
## #   ⁶​negpost_mean, ⁷​posbaseline_sd, ⁸​possession4_sd, ⁹​pospost_sd

Plots

#create id variable
ema_ac$id <- id_names
  
#reformat data: baseline
ema_ac_baseline <- ema_ac %>%
  dplyr::select(id, group, pos = posbaseline, neg = negbaseline)

#create time variable
ema_ac_baseline$time <- as.integer(1)
ema_ac_baseline <- relocate(ema_ac_baseline, time, .before=group)
#reformat data: session4
ema_ac_session4 <- ema_ac %>%
  dplyr::select(id, group, pos = possession4, neg = negsession4)

#create time variable
ema_ac_session4 $time <- as.integer(2)
ema_ac_session4  <- relocate(ema_ac_session4 , time, .before=group)
#reformat data: post
ema_ac_post <- ema_ac %>%
  dplyr::select(id, group, pos = pospost, neg = negpost)

#create time variable
ema_ac_post$time <- as.integer(3)
ema_ac_post<- relocate(ema_ac_post, time, .before=group)
#join predata and postdata
ema_ac2 <- merge(ema_ac_baseline, ema_ac_session4, all=TRUE)
ema_ac2 <- merge(ema_ac2, ema_ac_post, all=TRUE)

ema_ac2 $id <- factor(ema_ac2$id)
ema_ac2 $group <- factor(ema_ac2$group)
levels(ema_ac2 $group) <- c("Control", "Treatment")
ema_ac2 $time <- factor(ema_ac2$time)
levels(ema_ac2$time) <- c("baseline", "session4", "post")
#plot positive means
acpos_ggplot<- ggplot(data= ema_ac2, 
       aes(y=pos, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("AC of Positive Affect Scores Across Timepoints") + 
  labs(y = "autocorrelation", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
acpos_ggplot
## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Warning: Removed 1 row containing missing values (`geom_line()`).

#plot negative means
acneg_ggplot<- ggplot(data= ema_ac2, 
       aes(y=neg, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("AC of Negative Affect Scores Across Timepoints") + 
  labs(y = "autocorrelation", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
acneg_ggplot
## Warning: Removed 1 rows containing missing values (`geom_point()`).
## Warning: Removed 1 row containing missing values (`geom_line()`).

### Affect Instability (MSSD)

# initialize MSSD matrix
ema_mssd <- matrix(0, 32, 7)

colnames(ema_mssd) <- c("group", "posbaseline", "possession4", "pospost", "negbaseline", "negsession4", "negpost")
rownames(ema_mssd) <- id_names

for(i in 1:32){
  ema_mssd[i,1] <- id_group[i]
}
#fill positive mssd
## baseline
mssd <- mssd(betweenday_baseline$positive_mean, group = betweenday_baseline$id, lag = 1, na.rm = TRUE)

for(i in 1:32){
  ema_mssd[i,2] <- mssd[i]
}

## session4
mssd <- mssd(betweenday_session4$positive_mean, group = betweenday_session4$id, lag = 1, na.rm = TRUE)

for(i in 1:32){
  ema_mssd[i,3] <- mssd[i]
}

## post
mssd <- mssd(betweenday_post$positive_mean, group = betweenday_post$id, lag = 1, na.rm = TRUE)

for(i in 1:32){
  ema_mssd[i,4] <- mssd[i]
}
#fill negative mssd
## baseline
mssd <- mssd(betweenday_baseline$negative_mean, group = betweenday_baseline$id, lag = 1, na.rm = TRUE)

for(i in 1:32){
  ema_mssd[i,5] <- mssd[i]
}

## session 4
mssd <- mssd(betweenday_session4$negative_mean, group = betweenday_session4$id, lag = 1, na.rm = TRUE)

for(i in 1:32){
 ema_mssd[i,6] <- mssd[i]
}

## post
mssd <- mssd(betweenday_post$negative_mean, group = betweenday_post$id, lag = 1, na.rm = TRUE)

for(i in 1:32){
  ema_mssd[i,7] <- mssd[i]
}
ema_mssd <- as.data.frame(ema_mssd)
ema_mssd$posbaseline <- as.numeric(ema_mssd$posbaseline)
ema_mssd$possession4 <- as.numeric(ema_mssd$possession4)
ema_mssd$pospost <- as.numeric(ema_mssd$pospost)
ema_mssd$negbaseline <- as.numeric(ema_mssd$negbaseline)
ema_mssd$negsession4 <- as.numeric(ema_mssd$negsession4)
ema_mssd$negpost <- as.numeric(ema_mssd$negpost)
#group by group
ema_mssd%>%
  group_by(group) %>%
  summarise_at(vars(posbaseline, possession4, pospost, negbaseline, negsession4, negpost, ), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 2 × 13
##   group posbas…¹ posse…² pospo…³ negba…⁴ negse…⁵ negpo…⁶ posba…⁷ posse…⁸ pospo…⁹
##   <chr>    <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>
## 1 0         18.0    28.0    16.2    15.2    8.34   14.2     17.1    27.2    22.6
## 2 1         27.2    31.0    14.9    15.6    9.80    7.84    38.9    40.8    11.3
## # … with 3 more variables: negbaseline_sd <dbl>, negsession4_sd <dbl>,
## #   negpost_sd <dbl>, and abbreviated variable names ¹​posbaseline_mean,
## #   ²​possession4_mean, ³​pospost_mean, ⁴​negbaseline_mean, ⁵​negsession4_mean,
## #   ⁶​negpost_mean, ⁷​posbaseline_sd, ⁸​possession4_sd, ⁹​pospost_sd

Plots

#create id variable
ema_mssd$id <- id_names
  
#reformat data: baseline
ema_mssd_baseline <- ema_mssd %>%
  dplyr::select(id, group, pos = posbaseline, neg = negbaseline)

#create time variable
ema_mssd_baseline$time <- as.integer(1)
ema_mssd_baseline <- relocate(ema_mssd_baseline, time, .before=group)
#reformat data: session4
ema_mssd_session4 <- ema_mssd %>%
  dplyr::select(id, group, pos = possession4, neg = negsession4)

#create time variable
ema_mssd_session4 $time <- as.integer(2)
ema_mssd_session4  <- relocate(ema_mssd_session4 , time, .before=group)
#reformat data: post
ema_mssd_post <- ema_mssd %>%
  dplyr::select(id, group, pos = pospost, neg = negpost)

#create time variable
ema_mssd_post$time <- as.integer(3)
ema_mssd_post<- relocate(ema_mssd_post, time, .before=group)
#join predata and postdata
ema_mssd2 <- merge(ema_mssd_baseline, ema_mssd_session4, all=TRUE)
ema_mssd2 <- merge(ema_mssd2, ema_mssd_post, all=TRUE)

ema_mssd2 $id <- factor(ema_mssd2$id)
ema_mssd2 $group <- factor(ema_mssd2$group)
levels(ema_mssd2 $group) <- c("Control", "Treatment")
ema_mssd2 $time <- factor(ema_mssd2$time)
levels(ema_mssd2$time) <- c("baseline", "session4", "post")
#plot positive means
mssdpos_ggplot<- ggplot(data= ema_mssd2, 
       aes(y=pos, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("MSSD of Positive Affect Scores Across Timepoints") + 
  labs(y = "mssd", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
mssdpos_ggplot

#plot negative means
mssdneg_ggplot<- ggplot(data= ema_mssd2, 
       aes(y=neg, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("MSSD of Negative Affect Scores Across Timepoints") + 
  labs(y = "mssd", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
mssdneg_ggplot

Acute Increase (PAC)

# initialize PAC matrix
#baseline
ema_pac <- matrix(0, 32, 7)

colnames(ema_pac) <- c("group", "posbaseline", "possession4", "pospost", "negbaseline", "negsession4", "negpost")
rownames(ema_pac) <- id_names

for(i in 1:32){
  ema_pac[i,1] <- id_group[i]
}
#compute PAC for positive score
pac_vec <- c()
events <- c("Baseline", "Session 4", "Post")

for(n in 1:32){ #cycle through participant ids
  name <- id_names[n]
  
  for(e in 1:3){ #cycle through events
   
     index <- which(EMAdata$id == name & EMAdata$event == events[e])
      
    vec <- c()
    len <- length(index)
  
      for(i in 1:len){ #number of instances
        n <- EMAdata$positive_SD[index[i]]
        vec <- c(vec,n) #create vector of pos scores for each event
      
      }
      perc_90 <- quantile(vec, 0.9, na.rm=TRUE)
      ac <- ifelse(vec >= perc_90, 1, 0)
      pac <- sum(ac/(length(ac)-1),na.rm = TRUE)
    
      pac_vec <- c(pac_vec, pac)
    }
  }
pacpos <- matrix(pac_vec, nrow=32)

for(r in 1:32){
  ema_pac[r,2] <- pacpos[r,1]
  ema_pac[r,3] <- pacpos[r,2]
  ema_pac[r,4] <- pacpos[r,3]
}
#compute PAC for negative score
pac_vec <- c()
events <- c("Baseline", "Session 4", "Post")

for(n in 1:32){ #cycle through participant ids
  name <- id_names[n]
  
  for(e in 1:3){ #cycle through events
   
     index <- which(EMAdata$id == name & EMAdata$event == events[e])
      
    vec <- c()
    len <- length(index)
  
      for(i in 1:len){ #number of instances
        n <- EMAdata$negative_SD[index[i]]
        vec <- c(vec,n) #create vector of pos scores for each event
      
      }
      perc_90 <- quantile(vec, 0.9, na.rm=TRUE)
      ac <- ifelse(vec >= perc_90, 1, 0)
      pac <- sum(ac/(length(ac)-1),na.rm = TRUE)
    
      pac_vec <- c(pac_vec, pac)
    }
  }
pacneg <- matrix(pac_vec, nrow=32)

for(r in 1:32){
  ema_pac[r,5] <- pacneg[r,1]
  ema_pac[r,6] <- pacneg[r,2]
  ema_pac[r,7] <- pacneg[r,3]
}
ema_pac <- as.data.frame(ema_pac)
ema_pac$posbaseline <- as.numeric(ema_pac$posbaseline)
ema_pac$possession4 <- as.numeric(ema_pac$possession4)
ema_pac$pospost <- as.numeric(ema_pac$pospost)
ema_pac$negbaseline <- as.numeric(ema_pac$negbaseline)
ema_pac$negsession4 <- as.numeric(ema_pac$negsession4)
ema_pac$negpost <- as.numeric(ema_pac$negpost)
#group by group
ema_pac%>%
  group_by(group) %>%
  summarise_at(vars(posbaseline, possession4, pospost, negbaseline, negsession4, negpost, ), list(mean=mean, sd=sd), na.rm=TRUE)
## # A tibble: 2 × 13
##   group posbas…¹ posse…² pospo…³ negba…⁴ negse…⁵ negpo…⁶ posba…⁷ posse…⁸ pospo…⁹
##   <chr>    <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>
## 1 0        0.171   0.143   0.157   0.181   0.167   0.186  0.0591  0.0468  0.0400
## 2 1        0.172   0.160   0.181   0.168   0.155   0.164  0.0509  0.0312  0.0625
## # … with 3 more variables: negbaseline_sd <dbl>, negsession4_sd <dbl>,
## #   negpost_sd <dbl>, and abbreviated variable names ¹​posbaseline_mean,
## #   ²​possession4_mean, ³​pospost_mean, ⁴​negbaseline_mean, ⁵​negsession4_mean,
## #   ⁶​negpost_mean, ⁷​posbaseline_sd, ⁸​possession4_sd, ⁹​pospost_sd

Plots

#create id variable
ema_pac$id <- id_names

#reformat data: baseline
ema_pac_baseline <- ema_pac %>%
  dplyr::select(id, group, pos = posbaseline, neg = negbaseline)

#create time variable
ema_pac_baseline$time <- as.integer(1)
ema_pac_baseline <- relocate(ema_pac_baseline, time, .before=group)
#reformat data: session4
ema_pac_session4 <- ema_pac %>%
  dplyr::select(id, group, pos = possession4, neg = negsession4)

#create time variable
ema_pac_session4 $time <- as.integer(2)
ema_pac_session4  <- relocate(ema_pac_session4 , time, .before=group)
#reformat data: post
ema_pac_post <- ema_pac %>%
  dplyr::select(id, group, pos = pospost, neg = negpost)

#create time variable
ema_pac_post$time <- as.integer(3)
ema_pac_post<- relocate(ema_pac_post, time, .before=group)
#join predata and postdata
ema_pac2 <- merge(ema_pac_baseline, ema_pac_session4, all=TRUE)
ema_pac2 <- merge(ema_pac2, ema_pac_post, all=TRUE)

ema_pac2 $id <- factor(ema_pac2$id)
ema_pac2 $group <- factor(ema_pac2$group)
levels(ema_pac2 $group) <- c("Control", "Treatment")
ema_pac2 $time <- factor(ema_pac2$time)
levels(ema_pac2$time) <- c("baseline", "session4", "post")
#plot positive means
pacpos_ggplot<- ggplot(data= ema_pac2, 
       aes(y=pos, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("PAC of Positive Affect Scores Across Timepoints") + 
  labs(y = "mssd", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
pacpos_ggplot

#plot negative means
pacneg_ggplot<- ggplot(data= ema_pac2, 
       aes(y=neg, x=time, color=group, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("PAC of Negative Affect Scores Across Timepoints") + 
  labs(y = "mssd", x = "timepoint") +
  scale_color_manual(values = c("lightblue3","darkblue"))
 
pacneg_ggplot

Preliminary Analyses

mixed model w REML fixed effects: enrollment group, event, event*enrollmentgroup random effects: id

calculated in excel - kept running into issues doing it in R

filename <- "20230327_betweenday.csv"

betweenday <- read_csv(filename)
## Rows: 286 Columns: 16
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr   (1): date
## dbl  (12): id, event, day, instance, enrollmentgroup, completed, positive_sc...
## time  (3): time, time_participant, time_interval
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

x Events

lmm_emapos <- lmer(positive_mean ~ enrollmentgroup * event + (1|id), data=betweenday)

summary(lmm_emapos)
## Linear mixed model fit by REML ['lmerMod']
## Formula: positive_mean ~ enrollmentgroup * event + (1 | id)
##    Data: betweenday
## 
## REML criterion at convergence: 1729.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9201 -0.5531 -0.1101  0.5542  2.6795 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 24.04    4.903   
##  Residual             20.68    4.548   
## Number of obs: 283, groups:  id, 33
## 
## Fixed effects:
##                       Estimate Std. Error t value
## (Intercept)            21.9437     1.6064  13.660
## enrollmentgroup        -1.2525     2.2292  -0.562
## event                   0.6044     0.4806   1.258
## enrollmentgroup:event   0.1940     0.6660   0.291
## 
## Correlation of Fixed Effects:
##             (Intr) enrllm event 
## enrllmntgrp -0.721              
## event       -0.593  0.428       
## enrllmntgr:  0.428 -0.592 -0.722
anova(lmm_emapos) 
## Analysis of Variance Table
##                       npar Sum Sq Mean Sq F value
## enrollmentgroup          1  4.872   4.872  0.2356
## event                    1 92.984  92.984  4.4961
## enrollmentgroup:event    1  1.755   1.755  0.0849
eta_squared(lmm_emapos, partial=FALSE)
## Warning: Currently only supports partial eta squared for this class of objects.
## # Effect Size for ANOVA (Type III)
## 
## Parameter             | Eta2 (partial) |       95% CI
## -----------------------------------------------------
## enrollmentgroup       |       4.44e-03 | [0.00, 1.00]
## event                 |       6.32e-03 | [0.00, 1.00]
## enrollmentgroup:event |       3.41e-04 | [0.00, 1.00]
## 
## - One-sided CIs: upper bound fixed at [1.00].
lmm_emaneg <- lmer(negative_mean ~ enrollmentgroup * event + (1|id), data=betweenday)

summary(lmm_emaneg)
## Linear mixed model fit by REML ['lmerMod']
## Formula: negative_mean ~ enrollmentgroup * event + (1 | id)
##    Data: betweenday
## 
## REML criterion at convergence: 1361.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.4655 -0.6459 -0.0633  0.4641  4.3345 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 0.9872   0.9936  
##  Residual             6.5709   2.5634  
## Number of obs: 283, groups:  id, 33
## 
## Fixed effects:
##                       Estimate Std. Error t value
## (Intercept)           16.49769    0.63478  25.990
## enrollmentgroup        1.15568    0.87778   1.317
## event                 -0.03775    0.27048  -0.140
## enrollmentgroup:event -0.53124    0.37478  -1.417
## 
## Correlation of Fixed Effects:
##             (Intr) enrllm event 
## enrllmntgrp -0.723              
## event       -0.850  0.615       
## enrllmntgr:  0.613 -0.849 -0.722
anova(lmm_emaneg) 
## Analysis of Variance Table
##                       npar  Sum Sq Mean Sq F value
## enrollmentgroup          1  0.3255  0.3255  0.0495
## event                    1 18.5346 18.5346  2.8207
## enrollmentgroup:event    1 13.2022 13.2022  2.0092
eta_squared(lmm_emaneg, partial=FALSE)
## Warning: Currently only supports partial eta squared for this class of objects.
## # Effect Size for ANOVA (Type III)
## 
## Parameter             | Eta2 (partial) |       95% CI
## -----------------------------------------------------
## enrollmentgroup       |       7.90e-03 | [0.00, 1.00]
## event                 |       7.79e-05 | [0.00, 1.00]
## enrollmentgroup:event |       7.97e-03 | [0.00, 1.00]
## 
## - One-sided CIs: upper bound fixed at [1.00].

Plots

# Positive Affect plot
xlabels <- c("Baseline", "Session 4", "Post")

#create trendline
ci_panaspos <- confint(lmm_emapos, method="boot", nsim=10)
## Computing bootstrap confidence intervals ...
## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints
ci_panaspos <- tibble::rownames_to_column(data.frame(ci_panaspos), "Term") 
colnames(ci_panaspos)<- c("Term", "CI 2.5%", "CI 97.5%")
as.tibble(ci_panaspos)
## Warning: `as.tibble()` was deprecated in tibble 2.0.0.
## ℹ Please use `as_tibble()` instead.
## ℹ The signature and semantics have changed, see `?as_tibble`.
## # A tibble: 6 × 3
##   Term                  `CI 2.5%` `CI 97.5%`
##   <chr>                     <dbl>      <dbl>
## 1 .sig01                    4.41       5.81 
## 2 .sigma                    4.27       4.77 
## 3 (Intercept)              19.5       26.1  
## 4 enrollmentgroup          -8.90       2.07 
## 5 event                     0.254      1.23 
## 6 enrollmentgroup:event    -0.528      0.962
#Plot model - positive
panaspos_ggplot<- ggplot(data= betweenday, 
       aes(y=positive_mean, x=event, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Positive Affect Across Timepoints") + 
  labs(y = "positive affect score", x = "timepoint") +
  theme(legend.position = "none") + 
  scale_x_discrete(labels = xlabels) +
  geom_abline(intercept = fixef(lmm_emapos)[1], #Regression Line (RL).
              slope=fixef(lmm_emapos)[2], col="black", linetype = "dashed") +                 
  geom_abline(intercept = ci_panaspos$`CI 97.5%`[3], 
              slope=ci_panaspos$`CI 97.5%`[4], col="grey", linetype = "dotted") + #Upper Bound of RL
  geom_abline(intercept = ci_panaspos$`CI 2.5%`[3], 
              slope=ci_panaspos$`CI 2.5%`[4], col="grey", linetype= "dotted") #Lower Bound of RL
  
panaspos_ggplot
## Warning: Removed 3 rows containing missing values (`geom_point()`).
## Warning: Removed 3 rows containing missing values (`geom_line()`).

# Negative Affect plot
xlabels <- c("Baseline", "Session 4", "Post")

#create trendline
ci_panasneg <- confint(lmm_emaneg, method="boot", nsim=10)
## Computing bootstrap confidence intervals ...
## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints

## Warning in norm.inter(t, alpha): extreme order statistics used as endpoints
ci_panasneg <- tibble::rownames_to_column(data.frame(ci_panasneg), "Term") 
colnames(ci_panasneg)<- c("Term", "CI 2.5%", "CI 97.5%")
as.tibble(ci_panasneg)
## # A tibble: 6 × 3
##   Term                  `CI 2.5%` `CI 97.5%`
##   <chr>                     <dbl>      <dbl>
## 1 .sig01                   0.479      1.46  
## 2 .sigma                   2.50       2.74  
## 3 (Intercept)             15.1       17.4   
## 4 enrollmentgroup         -0.0586     2.75  
## 5 event                   -0.340      0.404 
## 6 enrollmentgroup:event   -1.11      -0.0456
#Plot model 
panasneg_ggplot<- ggplot(data= betweenday, 
       aes(y=negative_mean, x=event, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Negative Affect Across Timepoints") + 
  labs(y = "negative affect score", x = "timepoint") +
  theme(legend.position = "none") + 
  scale_x_discrete(labels = xlabels) +
  geom_abline(intercept = fixef(lmm_emaneg)[1], #Regression Line (RL).
              slope=fixef(lmm_emaneg)[2], col="black", linetype = "dashed") +                 
  geom_abline(intercept = ci_panasneg$`CI 97.5%`[3], 
              slope=ci_panasneg$`CI 97.5%`[4], col="grey", linetype = "dotted") + #Upper Bound of RL
  geom_abline(intercept = ci_panasneg$`CI 2.5%`[3], 
              slope=ci_panasneg$`CI 2.5%`[4], col="grey", linetype= "dotted") #Lower Bound of RL

panasneg_ggplot
## Warning: Removed 3 rows containing missing values (`geom_point()`).
## Warning: Removed 3 rows containing missing values (`geom_line()`).

x Days

betweenday$day <- as.factor(betweenday$day)
levels(betweenday$day) <- c("Day 1", "Day 2", "Day 3", "Day 4", "Day 5", "Day 6", "Day 7", "Day 8", "Day 9")
lmm_emapos <- lmer(positive_mean ~ enrollmentgroup * day + (1|id), data=betweenday)

anova(lmm_emapos) 
## Analysis of Variance Table
##                     npar  Sum Sq Mean Sq F value
## enrollmentgroup        1   4.865  4.8647  0.2369
## day                    8 236.891 29.6114  1.4420
## enrollmentgroup:day    8 185.041 23.1302  1.1264
summary(lmm_emapos)
## Linear mixed model fit by REML ['lmerMod']
## Formula: positive_mean ~ enrollmentgroup * day + (1 | id)
##    Data: betweenday
## 
## REML criterion at convergence: 1676
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -2.53516 -0.57192 -0.01372  0.53747  2.70835 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 23.92    4.891   
##  Residual             20.53    4.531   
## Number of obs: 283, groups:  id, 33
## 
## Fixed effects:
##                          Estimate Std. Error t value
## (Intercept)               22.8416     1.6514  13.831
## enrollmentgroup            0.7329     2.3113   0.317
## dayDay 2                  -0.8544     1.6393  -0.521
## dayDay 3                  -0.6758     1.6722  -0.404
## dayDay 4                   1.2023     1.6393   0.733
## dayDay 5                   0.9223     1.6393   0.563
## dayDay 6                   0.2023     1.6393   0.123
## dayDay 7                  -0.2288     1.6393  -0.140
## dayDay 8                   0.9823     1.6393   0.599
## dayDay 9                   1.2245     1.6393   0.747
## enrollmentgroup:dayDay 2  -1.1662     2.2590  -0.516
## enrollmentgroup:dayDay 3  -1.4060     2.3023  -0.611
## enrollmentgroup:dayDay 4  -3.3680     2.2590  -1.491
## enrollmentgroup:dayDay 5  -5.0058     2.2776  -2.198
## enrollmentgroup:dayDay 6  -2.3069     2.2776  -1.013
## enrollmentgroup:dayDay 7  -0.3644     2.2785  -0.160
## enrollmentgroup:dayDay 8   0.4162     2.2785   0.183
## enrollmentgroup:dayDay 9  -1.1573     2.2785  -0.508
## 
## Correlation matrix not shown by default, as p = 18 > 12.
## Use print(x, correlation=TRUE)  or
##     vcov(x)        if you need it
eta_squared(lmm_emapos, partial = FALSE)
## Warning: Currently only supports partial eta squared for this class of objects.
## # Effect Size for ANOVA (Type III)
## 
## Parameter           | Eta2 (partial) |       95% CI
## ---------------------------------------------------
## enrollmentgroup     |       7.43e-03 | [0.00, 1.00]
## day                 |           0.02 | [0.00, 1.00]
## enrollmentgroup:day |           0.04 | [0.00, 1.00]
## 
## - One-sided CIs: upper bound fixed at [1.00].
lmm_emaneg <- lmer(negative_mean ~ enrollmentgroup * day + (1|id), data=betweenday)

anova(lmm_emaneg) 
## Analysis of Variance Table
##                     npar Sum Sq Mean Sq F value
## enrollmentgroup        1  0.330  0.3295  0.0503
## day                    8 75.042  9.3803  1.4310
## enrollmentgroup:day    8 54.495  6.8118  1.0392
summary(lmm_emaneg)
## Linear mixed model fit by REML ['lmerMod']
## Formula: negative_mean ~ enrollmentgroup * day + (1 | id)
##    Data: betweenday
## 
## REML criterion at convergence: 1325
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.3481 -0.6159 -0.0443  0.4982  4.5235 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 0.9708   0.9853  
##  Residual             6.5552   2.5603  
## Number of obs: 283, groups:  id, 33
## 
## Fixed effects:
##                          Estimate Std. Error t value
## (Intercept)               17.0253     0.6699  25.416
## enrollmentgroup           -0.1939     0.9442  -0.205
## dayDay 2                  -0.4528     0.9131  -0.496
## dayDay 3                  -0.6556     0.9312  -0.704
## dayDay 4                  -0.7484     0.9131  -0.820
## dayDay 5                  -1.2128     0.9131  -1.328
## dayDay 6                  -1.2306     0.9131  -1.348
## dayDay 7                  -0.3839     0.9131  -0.420
## dayDay 8                  -0.7661     0.9131  -0.839
## dayDay 9                  -0.1261     0.9131  -0.138
## enrollmentgroup:dayDay 2   0.3254     1.2669   0.257
## enrollmentgroup:dayDay 3   2.2505     1.2902   1.744
## enrollmentgroup:dayDay 4   0.5052     1.2669   0.399
## enrollmentgroup:dayDay 5   0.6789     1.2770   0.532
## enrollmentgroup:dayDay 6  -0.2686     1.2770  -0.210
## enrollmentgroup:dayDay 7   0.3049     1.2772   0.239
## enrollmentgroup:dayDay 8   0.3694     1.2772   0.289
## enrollmentgroup:dayDay 9  -1.3633     1.2772  -1.067
## 
## Correlation matrix not shown by default, as p = 18 > 12.
## Use print(x, correlation=TRUE)  or
##     vcov(x)        if you need it
eta_squared(lmm_emaneg, partial = FALSE)
## Warning: Currently only supports partial eta squared for this class of objects.
## # Effect Size for ANOVA (Type III)
## 
## Parameter           | Eta2 (partial) |       95% CI
## ---------------------------------------------------
## enrollmentgroup     |       2.07e-03 | [0.00, 1.00]
## day                 |           0.01 | [0.00, 1.00]
## enrollmentgroup:day |           0.03 | [0.00, 1.00]
## 
## - One-sided CIs: upper bound fixed at [1.00].

Plots

# Positive Affect plot
xlabels <- c("Day 1", "Day 2", "Day 3", "Day 4", "Day 5", "Day 6", "Day 7", "Day 8", "Day 9")

#Plot model - positive
panaspos_ggplot<- ggplot(data= betweenday, 
       aes(y=positive_mean, x=day, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Positive Affect Across Days") + 
  labs(y = "positive affect score", x = "day") +
  theme(legend.position = "none") + 
  scale_x_discrete(labels = xlabels) 
 
panaspos_ggplot
## Warning: Removed 3 rows containing missing values (`geom_point()`).
## Warning: Removed 3 rows containing missing values (`geom_line()`).

# Negative Affect plot
xlabels <- c("Day 1", "Day 2", "Day 3", "Day 4", "Day 5", "Day 6", "Day 7", "Day 8", "Day 9")

#Plot model - positive
panaspos_ggplot<- ggplot(data= betweenday, 
       aes(y=negative_mean, x=day, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Negative Affect Across Days") + 
  labs(y = "negative affect score", x = "day") +
  theme(legend.position = "none") + 
  scale_x_discrete(labels = xlabels)

panaspos_ggplot
## Warning: Removed 3 rows containing missing values (`geom_point()`).
## Warning: Removed 3 rows containing missing values (`geom_line()`).

x Instances

EMAdata$instance <- as.factor(EMAdata$instance)
levels(EMAdata$instance) <- c("1","2","3","4","5","6","7","8","9","10","11","12","13","14","15")
lmm_emapos <- lmer(positive_score ~ enrollmentgroup * instance + (1|id), data= EMAdata)

anova(lmm_emapos) 
## Analysis of Variance Table
##                          npar  Sum Sq Mean Sq F value
## enrollmentgroup             1    0.89   0.887  0.0268
## instance                   14 2986.95 213.353  6.4483
## enrollmentgroup:instance   14  473.82  33.845  1.0229
summary(lmm_emapos)
## Linear mixed model fit by REML ['lmerMod']
## Formula: positive_score ~ enrollmentgroup * instance + (1 | id)
##    Data: EMAdata
## 
## REML criterion at convergence: 7896.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.6314 -0.5667 -0.0795  0.5492  4.7119 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 27.08    5.204   
##  Residual             33.09    5.752   
## Number of obs: 1242, groups:  id, 33
## 
## Fixed effects:
##                                                Estimate Std. Error t value
## (Intercept)                                    24.96158    1.62308  15.379
## enrollmentgroupImmediate Enrollment            -3.06792    2.22861  -1.377
## instance2                                       0.03676    1.29416   0.028
## instance3                                      -1.31717    1.32381  -0.995
## instance4                                      -2.50545    1.29935  -1.928
## instance5                                      -5.40576    1.30753  -4.134
## instance6                                      -2.42108    1.30664  -1.853
## instance7                                       0.46633    1.31510   0.355
## instance8                                      -1.14320    1.29980  -0.880
## instance9                                      -2.95746    1.28561  -2.300
## instance10                                     -3.91047    1.34222  -2.913
## instance11                                     -1.43807    1.32299  -1.087
## instance12                                     -1.27751    1.31507  -0.971
## instance13                                     -1.30025    1.29934  -1.001
## instance14                                     -1.55438    1.31513  -1.182
## instance15                                     -5.22615    1.31547  -3.973
## enrollmentgroupImmediate Enrollment:instance2   2.71039    1.78096   1.522
## enrollmentgroupImmediate Enrollment:instance3   3.26754    1.81916   1.796
## enrollmentgroupImmediate Enrollment:instance4   4.25859    1.80876   2.354
## enrollmentgroupImmediate Enrollment:instance5   4.00111    1.81513   2.204
## enrollmentgroupImmediate Enrollment:instance6   2.63457    1.77656   1.483
## enrollmentgroupImmediate Enrollment:instance7   2.57019    1.79726   1.430
## enrollmentgroupImmediate Enrollment:instance8   1.39434    1.79168   0.778
## enrollmentgroupImmediate Enrollment:instance9   3.87343    1.80248   2.149
## enrollmentgroupImmediate Enrollment:instance10  1.74580    1.82771   0.955
## enrollmentgroupImmediate Enrollment:instance11  0.99683    1.80811   0.551
## enrollmentgroupImmediate Enrollment:instance12  4.11750    1.81197   2.272
## enrollmentgroupImmediate Enrollment:instance13  3.99321    1.78039   2.243
## enrollmentgroupImmediate Enrollment:instance14  2.85621    1.81530   1.573
## enrollmentgroupImmediate Enrollment:instance15  3.76380    1.82505   2.062
## 
## Correlation matrix not shown by default, as p = 30 > 12.
## Use print(x, correlation=TRUE)  or
##     vcov(x)        if you need it
lme_emaneg <- lmer(negative_score ~ enrollmentgroup * instance + (1|id), data= EMAdata)

anova(lme_emaneg) 
## Analysis of Variance Table
##                          npar Sum Sq Mean Sq F value
## enrollmentgroup             1   1.20   1.197  0.0942
## instance                   14 601.54  42.967  3.3811
## enrollmentgroup:instance   14  90.62   6.473  0.5093
summary(lme_emaneg)
## Linear mixed model fit by REML ['lmerMod']
## Formula: negative_score ~ enrollmentgroup * instance + (1 | id)
##    Data: EMAdata
## 
## REML criterion at convergence: 6730.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.6836 -0.5522 -0.1338  0.3839  6.0894 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.447   2.906   
##  Residual             12.708   3.565   
## Number of obs: 1242, groups:  id, 33
## 
## Fixed effects:
##                                                Estimate Std. Error t value
## (Intercept)                                    16.96455    0.94258  17.998
## enrollmentgroupImmediate Enrollment            -0.27901    1.29204  -0.216
## instance2                                       0.37199    0.80162   0.464
## instance3                                      -0.50098    0.82041  -0.611
## instance4                                      -0.73212    0.80526  -0.909
## instance5                                      -1.40252    0.81033  -1.731
## instance6                                      -0.63540    0.80978  -0.785
## instance7                                      -0.41365    0.81502  -0.508
## instance8                                      -0.45366    0.80554  -0.563
## instance9                                      -1.01297    0.79675  -1.271
## instance10                                     -1.79199    0.83182  -2.154
## instance11                                      0.01560    0.81991   0.019
## instance12                                     -0.33090    0.81500  -0.406
## instance13                                     -1.66764    0.80525  -2.071
## instance14                                     -1.32012    0.81504  -1.620
## instance15                                     -1.86610    0.81525  -2.289
## enrollmentgroupImmediate Enrollment:instance2   0.04248    1.10341   0.039
## enrollmentgroupImmediate Enrollment:instance3  -1.18564    1.12738  -1.052
## enrollmentgroupImmediate Enrollment:instance4  -0.58221    1.12092  -0.519
## enrollmentgroupImmediate Enrollment:instance5  -0.53164    1.12486  -0.473
## enrollmentgroupImmediate Enrollment:instance6  -0.25795    1.10099  -0.234
## enrollmentgroupImmediate Enrollment:instance7   0.47638    1.11381   0.428
## enrollmentgroupImmediate Enrollment:instance8  -0.36636    1.11034  -0.330
## enrollmentgroupImmediate Enrollment:instance9   0.60958    1.11705   0.546
## enrollmentgroupImmediate Enrollment:instance10  0.35344    1.13267   0.312
## enrollmentgroupImmediate Enrollment:instance11  0.02329    1.12053   0.021
## enrollmentgroupImmediate Enrollment:instance12  0.49568    1.12293   0.441
## enrollmentgroupImmediate Enrollment:instance13  0.98644    1.10336   0.894
## enrollmentgroupImmediate Enrollment:instance14 -0.35677    1.12497  -0.317
## enrollmentgroupImmediate Enrollment:instance15 -0.36472    1.13103  -0.322
## 
## Correlation matrix not shown by default, as p = 30 > 12.
## Use print(x, correlation=TRUE)  or
##     vcov(x)        if you need it

Plots

# Positive Affect plot @ Baseline
EMA_baseline <- filter(EMAdata, event == "Baseline")
  
#Plot model - positive
instance_pos<- ggplot(data= EMA_baseline, 
       aes(y=positive_score, x=instance, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Positive Affect Across Instance @ Baseline") + 
  labs(y = "positive affect score", x = "instance") +
  theme(legend.position = "none") 

instance_pos
## Warning: Removed 82 rows containing missing values (`geom_point()`).
## Warning: Removed 24 rows containing missing values (`geom_line()`).

# Positive Affect plot @ Session 4
EMA_session4 <- filter(EMAdata, event == "Session 4")
  
#Plot model - positive
instance_pos<- ggplot(data= EMA_session4, 
       aes(y=positive_score, x=instance, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Positive Affect Across Instance @ Session4") + 
  labs(y = "positive affect score", x = "instance") +
  theme(legend.position = "none") 

instance_pos
## Warning: Removed 44 rows containing missing values (`geom_point()`).
## Warning: Removed 16 rows containing missing values (`geom_line()`).

# Positive Affect plot @ Post
EMA_post <- filter(EMAdata, event == "Post")
  
#Plot model - positive
instance_pos<- ggplot(data= EMA_post, 
       aes(y=positive_score, x=instance, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Positive Affect Across Instance @ Post") + 
  labs(y = "positive affect score", x = "instance") +
  theme(legend.position = "none") 

instance_pos
## Warning: Removed 59 rows containing missing values (`geom_point()`).
## Warning: Removed 17 rows containing missing values (`geom_line()`).

# Negative Affect plot @ Baseline

instance_neg<- ggplot(data= EMA_baseline, 
       aes(y=negative_score, x=instance, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Negative Affect Across Instance  @ Baseline") + 
  labs(y = "negative affect score", x = "instance") +
  theme(legend.position = "none") 

instance_neg
## Warning: Removed 82 rows containing missing values (`geom_point()`).
## Warning: Removed 24 rows containing missing values (`geom_line()`).

# Negative Affect plot @ Session 4

instance_neg<- ggplot(data= EMA_session4, 
       aes(y=negative_score, x=instance, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Negative Affect Across Instance  @ Session 4") + 
  labs(y = "negative affect score", x = "instance") +
  theme(legend.position = "none") 

instance_neg
## Warning: Removed 44 rows containing missing values (`geom_point()`).
## Warning: Removed 16 rows containing missing values (`geom_line()`).

# Negative Affect plot @ Post

instance_neg<- ggplot(data= EMA_post, 
       aes(y=negative_score, x=instance, color=id, group=id)) + 
  geom_point() + geom_line() +
  ggtitle("Negative Affect Across Instance  @ Post") + 
  labs(y = "negative affect score", x = "instance") +
  theme(legend.position = "none") 

instance_neg
## Warning: Removed 59 rows containing missing values (`geom_point()`).
## Warning: Removed 17 rows containing missing values (`geom_line()`).

Multilevel Model

#combine dataset 
mlm <- ema_mean2[,1:3]

mlm$pos_mean <- ema_mean2[,4]
mlm$neg_mean <- ema_mean2[,5]

mlm$pos_sd <- ema_sd2[,4]
mlm$neg_sd <- ema_sd2[,5]

mlm$pos_ac <- ema_ac2[,4]
mlm$neg_ac <- ema_ac2[,5]

mlm$pos_mssd <- ema_mssd2[,4]
mlm$neg_mssd <- ema_mssd2[,5]

mlm$pos_pac <- ema_pac2[,4]
mlm$neg_pac <- ema_pac2[,5]
outcome <- read_csv("20230325_phase2data.csv")
## Rows: 108 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): event
## dbl (8): id, group, ryff, pss, bdi, bai, fmqi, cerq
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
outcome <- filter(outcome, group != 2)
outcome <- filter(outcome, id != "12046")
outcome <- filter(outcome, event == "baseline") #baseline only

outcome <- outcome[rep(seq_len(nrow(outcome)), each = 3), ]

class(outcome)
## [1] "tbl_df"     "tbl"        "data.frame"
outc <- as.data.frame(outcome[,4:8])
class(outc)
## [1] "data.frame"
mlm$ryff <- (outc[,1])
mlm$pss <- outc[,2]
mlm$bdi <- outc[,3]
mlm$bai <- outc[,4]
mlm$mi <- outc[,5]

view(mlm)
write.csv(mlm, "20230411_mlm.csv")
mlm <- read.csv("MLM_USETHIS.csv")
library(lmerTest) #tests on lmer objects
## 
## Attaching package: 'lmerTest'
## The following object is masked from 'package:lme4':
## 
##     lmer
## The following object is masked from 'package:stats':
## 
##     step
library(performance) # ICC calculations
library(texreg)
## Version:  1.38.6
## Date:     2022-04-06
## Author:   Philip Leifeld (University of Essex)
## 
## Consider submitting praise using the praise or praise_interactive functions.
## Please cite the JSS article in your publications -- see citation("texreg").
## 
## Attaching package: 'texreg'
## The following object is masked from 'package:tidyr':
## 
##     extract
#library(interactions) # plotting interactions

Null Model

#Random-Intercept Only
nullmodel_pos <- lmer(positive_score ~ 1 + (1|id), data = mlm)
performance::icc(nullmodel_pos)
## # Intraclass Correlation Coefficient
## 
##     Adjusted ICC: 0.433
##   Unadjusted ICC: 0.433
#Random-Intercept Only
nullmodel_neg <- lmer(negative_score ~ 1 + (1|id), data = mlm)
performance::icc(nullmodel_neg)
## # Intraclass Correlation Coefficient
## 
##     Adjusted ICC: 0.393
##   Unadjusted ICC: 0.393

Level-1 (Instance)

Fixed Effects

#L1 Fixed Effects

l1model_pos <- lmer(positive_score ~ 1 + instance + (1|id), data = mlm)
summary(l1model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 8040.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.6373 -0.5897 -0.0500  0.6045  4.2817 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 26.72    5.169   
##  Residual             35.09    5.924   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   23.41100    0.98062   38.17234  23.874   <2e-16 ***
## instance      -0.08331    0.03917 1207.66345  -2.127   0.0336 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr)
## instance -0.316
anova(nullmodel_pos, l1model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## nullmodel_pos: positive_score ~ 1 + (1 | id)
## l1model_pos: positive_score ~ 1 + instance + (1 | id)
##               npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)  
## nullmodel_pos    3 8047.9 8063.2 -4020.9   8041.9                       
## l1model_pos      4 8045.3 8065.8 -4018.7   8037.3 4.5243  1    0.03342 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#L1 Fixed Effects

l1model_neg <- lmer(negative_score ~ 1 + instance + (1|id), data = mlm)
summary(l1model_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6798
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.0042 -0.5198 -0.1522  0.3655  6.3319 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.393   2.897   
##  Residual             12.913   3.594   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   16.55811    0.55573   39.68587   29.80  < 2e-16 ***
## instance      -0.06867    0.02376 1207.96661   -2.89  0.00392 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr)
## instance -0.338
anova(nullmodel_neg, l1model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## nullmodel_neg: negative_score ~ 1 + (1 | id)
## l1model_neg: negative_score ~ 1 + instance + (1 | id)
##               npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)   
## nullmodel_neg    3 6807.2 6822.5 -3400.6   6801.2                        
## l1model_neg      4 6800.8 6821.3 -3396.4   6792.8 8.3338  1   0.003891 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Random Slopes

l1random_pos <- lmer(positive_score ~ 1 + instance + (1|id) + 
                       (0 + instance|id), data= mlm)
summary(l1random_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + (1 | id) + (0 + instance | id)
##    Data: mlm
## 
## REML criterion at convergence: 8023.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.5319 -0.5766 -0.0284  0.5911  4.3979 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 22.8910  4.7845  
##  id.1     instance     0.0633  0.2516  
##  Residual             33.9445  5.8262  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##             Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept) 23.36757    0.91618 34.23047  25.505   <2e-16 ***
## instance    -0.07392    0.05921 32.96449  -1.248    0.221    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr)
## instance -0.221
anova(l1random_pos, l1model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l1model_pos: positive_score ~ 1 + instance + (1 | id)
## l1random_pos: positive_score ~ 1 + instance + (1 | id) + (0 + instance | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l1model_pos     4 8045.3 8065.8 -4018.7   8037.3                         
## l1random_pos    5 8031.4 8057.1 -4010.7   8021.4 15.894  1  6.699e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
l1random_neg <- lmer(negative_score ~ 1 + instance + (1|id) + 
                       (0 + instance|id), data= mlm)
summary(l1random_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + (1 | id) + (0 + instance | id)
##    Data: mlm
## 
## REML criterion at convergence: 6783.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5960 -0.5264 -0.1387  0.3486  6.0853 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.8866  2.9810  
##  id.1     instance     0.0227  0.1507  
##  Residual             12.4519  3.5287  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##             Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept) 16.52490    0.56852 33.81751  29.067   <2e-16 ***
## instance    -0.06135    0.03563 30.56734  -1.722   0.0952 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr)
## instance -0.217
anova(l1random_neg, l1model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l1model_neg: negative_score ~ 1 + instance + (1 | id)
## l1random_neg: negative_score ~ 1 + instance + (1 | id) + (0 + instance | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l1model_neg     4 6800.8 6821.3 -3396.4   6792.8                         
## l1random_neg    5 6789.1 6814.7 -3389.6   6779.1 13.731  1  0.0002109 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Level-2 (Day)

Fixed Effects

l2model_pos <- lmer(positive_score ~ 1 + instance + day +
                      (1|id) + (0 + instance|id), data= mlm)
summary(l2model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + (1 | id) + (0 + instance |  
##     id)
##    Data: mlm
## 
## REML criterion at convergence: 7994.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.5174 -0.5429 -0.0335  0.5527  4.6759 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 22.49770 4.7432  
##  id.1     instance     0.06287 0.2507  
##  Residual             33.19207 5.7613  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)   21.4594     0.9755   45.4452  21.999  < 2e-16 ***
## instance      -0.6651     0.1252  489.7917  -5.313 1.64e-07 ***
## day            3.3056     0.6178 1184.3983   5.350 1.05e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc
## instance  0.226       
## day      -0.366 -0.883
anova(l1random_pos, l2model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l1random_pos: positive_score ~ 1 + instance + (1 | id) + (0 + instance | id)
## l2model_pos: positive_score ~ 1 + instance + day + (1 | id) + (0 + instance | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l1random_pos    5 8031.4 8057.1 -4010.7   8021.4                         
## l2model_pos     6 8005.1 8035.8 -3996.5   7993.1 28.351  1  1.012e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
l2model_neg <- lmer(negative_score ~ 1 + instance + day +
                      (1|id) + (0 + instance|id), data= mlm)
summary(l2model_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + (1 | id) + (0 + instance |  
##     id)
##    Data: mlm
## 
## REML criterion at convergence: 6757.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5588 -0.5383 -0.1254  0.3518  5.9629 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.90385 2.9839  
##  id.1     instance     0.02291 0.1514  
##  Residual             12.19419 3.4920  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   15.42320    0.60790   44.10511  25.371  < 2e-16 ***
## instance      -0.40255    0.07584  475.56427  -5.308 1.71e-07 ***
## day            1.90816    0.37450 1172.25636   5.095 4.06e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc
## instance  0.220       
## day      -0.356 -0.883
anova(l1random_neg, l2model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l1random_neg: negative_score ~ 1 + instance + (1 | id) + (0 + instance | id)
## l2model_neg: negative_score ~ 1 + instance + day + (1 | id) + (0 + instance | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l1random_neg    5 6789.1 6814.7 -3389.6   6779.1                         
## l2model_neg     6 6765.4 6796.1 -3376.7   6753.4 25.708  1  3.972e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Random Slopes

l2random_pos <- lmer(positive_score ~ 1 + instance + day +
                      (1|id) + (0 + instance|id) + (0 + day|id), 
                     data= mlm)
## boundary (singular) fit: see help('isSingular')
summary(l2random_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + (1 | id) + (0 + instance |  
##     id) + (0 + day | id)
##    Data: mlm
## 
## REML criterion at convergence: 7994.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.5174 -0.5429 -0.0335  0.5527  4.6759 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 22.49152 4.7425  
##  id.1     instance     0.06286 0.2507  
##  id.2     day          0.00000 0.0000  
##  Residual             33.19235 5.7613  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)   21.4594     0.9754   45.4656  22.001  < 2e-16 ***
## instance      -0.6651     0.1252  489.8648  -5.313 1.64e-07 ***
## day            3.3056     0.6179 1184.3939   5.350 1.05e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc
## instance  0.227       
## day      -0.366 -0.883
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
anova(l2model_pos, l2random_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l2model_pos: positive_score ~ 1 + instance + day + (1 | id) + (0 + instance | id)
## l2random_pos: positive_score ~ 1 + instance + day + (1 | id) + (0 + instance | id) + (0 + day | id)
##              npar    AIC    BIC  logLik deviance Chisq Df Pr(>Chisq)
## l2model_pos     6 8005.1 8035.8 -3996.5   7993.1                    
## l2random_pos    7 8007.1 8043.0 -3996.5   7993.1     0  1          1
l2random_neg <- lmer(negative_score ~ 1 + instance + day +
                      (1|id) + (0 + instance|id) + (0 + day|id), 
                     data= mlm)

summary(l2random_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + (1 | id) + (0 + instance |  
##     id) + (0 + day | id)
##    Data: mlm
## 
## REML criterion at convergence: 6752.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.4727 -0.5293 -0.1162  0.3489  6.0696 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  8.724935 2.95380 
##  id.1     instance     0.005796 0.07613 
##  id.2     day          0.606159 0.77856 
##  Residual             12.088054 3.47679 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)  15.40876    0.60331  39.90184  25.541  < 2e-16 ***
## instance     -0.40666    0.07188 323.87748  -5.658 3.38e-08 ***
## day           1.92933    0.39827 383.75487   4.844 1.85e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc
## instance  0.235       
## day      -0.338 -0.870
anova(l2model_neg, l2random_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l2model_neg: negative_score ~ 1 + instance + day + (1 | id) + (0 + instance | id)
## l2random_neg: negative_score ~ 1 + instance + day + (1 | id) + (0 + instance | id) + (0 + day | id)
##              npar    AIC    BIC  logLik deviance Chisq Df Pr(>Chisq)  
## l2model_neg     6 6765.4 6796.1 -3376.7   6753.4                      
## l2random_neg    7 6761.9 6797.8 -3374.0   6747.9  5.47  1    0.01935 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Level-3 (Event)

Fixed Effects

l3model_pos <- lmer(positive_score ~ 1 + instance + day + event +
                      (1|id) + (0 + instance|id), data= mlm)

summary(l3model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + (1 | id) + (0 +  
##     instance | id)
##    Data: mlm
## 
## REML criterion at convergence: 7950.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3328 -0.6053 -0.0121  0.5105  4.5273 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 22.50731 4.7442  
##  id.1     instance     0.06286 0.2507  
##  Residual             31.98979 5.6560  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)   18.7438     1.0508   61.3227  17.838  < 2e-16 ***
## instance      -0.6497     0.1232  476.7678  -5.273 2.04e-07 ***
## day            3.2475     0.6066 1182.9940   5.353 1.04e-07 ***
## event          1.3569     0.2007 1179.7832   6.763 2.13e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day   
## instance  0.199              
## day      -0.328 -0.881       
## event    -0.382  0.018 -0.014
anova(l2model_pos, l3model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l2model_pos: positive_score ~ 1 + instance + day + (1 | id) + (0 + instance | id)
## l3model_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id)
##             npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l2model_pos    6 8005.1 8035.8 -3996.5   7993.1                         
## l3model_pos    7 7962.1 7998.0 -3974.1   7948.1 44.973  1  1.998e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
l3model_neg <- lmer(negative_score ~ 1 + instance + day + event +
                      (1|id) + (0 + instance|id) + (0 + day|id), 
                     data= mlm)

summary(l3model_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + (1 | id) + (0 +  
##     instance | id) + (0 + day | id)
##    Data: mlm
## 
## REML criterion at convergence: 6751.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.5400 -0.5506 -0.1118  0.3503  6.0143 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  8.731818 2.95496 
##  id.1     instance     0.005723 0.07565 
##  id.2     day          0.607983 0.77973 
##  Residual             12.065446 3.47354 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   15.85093    0.65180   54.32332  24.319  < 2e-16 ***
## instance      -0.40906    0.07181  323.80669  -5.697 2.75e-08 ***
## day            1.93852    0.39805  383.52098   4.870 1.63e-06 ***
## event         -0.22096    0.12323 1150.14318  -1.793   0.0732 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day   
## instance  0.210              
## day      -0.307 -0.870       
## event    -0.378  0.019 -0.013
anova(l2random_neg, l3model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l2random_neg: negative_score ~ 1 + instance + day + (1 | id) + (0 + instance | id) + (0 + day | id)
## l3model_neg: negative_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + day | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)  
## l2random_neg    7 6761.9 6797.8 -3374.0   6747.9                       
## l3model_neg     8 6760.7 6801.7 -3372.4   6744.7 3.2159  1    0.07292 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Random Effects

l3random_pos <- lmer(positive_score ~ 1 + instance + day + event +
                      (1|id) + (0 + instance|id) + (0 + event|id), 
                     data= mlm)

summary(l3random_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + (1 | id) + (0 +  
##     instance | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 7792.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.9869 -0.5841 -0.0094  0.5641  4.0528 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 39.99764 6.3244  
##  id.1     instance     0.06029 0.2455  
##  id.2     event        7.76956 2.7874  
##  Residual             26.09718 5.1085  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)   18.9336     1.2602   36.0397  15.025  < 2e-16 ***
## instance      -0.6696     0.1128  421.0196  -5.934 6.17e-09 ***
## day            3.3404     0.5494 1137.1313   6.080 1.64e-09 ***
## event          1.2054     0.5288   31.1688   2.280   0.0296 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day   
## instance  0.148              
## day      -0.246 -0.871       
## event    -0.108  0.007 -0.007
anova(l3model_pos, l3random_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3model_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id)
## l3random_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l3model_pos     7 7962.1 7998.0 -3974.1   7948.1                         
## l3random_pos    8 7807.6 7848.6 -3895.8   7791.6 156.51  1  < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
l3random_neg <- lmer(negative_score ~ 1 + instance + day + event +
                      (1|id) + (0 + instance|id) + (0 + day|id) + 
                       (0 + event|id), data= mlm)

summary(l3random_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + (1 | id) + (0 +  
##     instance | id) + (0 + day | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 6692.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7614 -0.5512 -0.1165  0.3426  5.8309 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  8.887714 2.98123 
##  id.1     instance     0.008473 0.09205 
##  id.2     day          0.852399 0.92325 
##  id.3     event        1.430644 1.19610 
##  Residual             10.951782 3.30935 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)  16.01054    0.64699  38.08214  24.746  < 2e-16 ***
## instance     -0.41459    0.06935 253.77099  -5.978 7.63e-09 ***
## day           1.97995    0.39254 291.02557   5.044 8.04e-07 ***
## event        -0.33384    0.24423  30.00279  -1.367    0.182    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day   
## instance  0.199              
## day      -0.285 -0.834       
## event    -0.185  0.012 -0.010
anova(l3model_neg, l3random_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3model_neg: negative_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + day | id)
## l3random_neg: negative_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l3model_neg     8 6760.7 6801.7 -3372.4   6744.7                         
## l3random_neg    9 6705.7 6751.8 -3343.9   6687.7 56.989  1  4.383e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Level-4 (Participant)

Fixed Effects

l4model_pos <- lmer(positive_score ~ 1 + instance + day + event + group +
                      (1|id) + (0 + instance|id) + (0 + event|id), 
                     data= mlm)

summary(l4model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + group + (1 | id) +  
##     (0 + instance | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 7788.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.9870 -0.5836 -0.0093  0.5637  4.0520 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 41.4187  6.4357  
##  id.1     instance     0.0607  0.2464  
##  id.2     event        7.8084  2.7944  
##  Residual             26.0912  5.1080  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)   19.2834     1.8136   30.8734  10.633 7.63e-12 ***
## instance      -0.6696     0.1129  418.3341  -5.932 6.29e-09 ***
## day            3.3406     0.5494 1137.6212   6.081 1.63e-09 ***
## event          1.2056     0.5300   31.1639   2.275   0.0299 *  
## group         -0.6591     2.4440   28.5776  -0.270   0.7894    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event 
## instance  0.101                     
## day      -0.168 -0.871              
## event    -0.068  0.007 -0.007       
## group    -0.710  0.002 -0.005 -0.009
anova(l3random_pos, l4model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + event | id)
## l4model_pos: positive_score ~ 1 + instance + day + event + group + (1 | id) + (0 + instance | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)
## l3random_pos    8 7807.6 7848.6 -3895.8   7791.6                     
## l4model_pos     9 7809.5 7855.6 -3895.8   7791.5 0.0797  1     0.7777
l4model_neg <- lmer(negative_score ~ 1 + instance + day + event + group +
                      (1|id) + (0 + instance|id) + (0 + day|id) + 
                       (0 + event|id), data= mlm)

summary(l4model_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + group + (1 | id) +  
##     (0 + instance | id) + (0 + day | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 6690.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7640 -0.5471 -0.1162  0.3393  5.8382 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  9.173949 3.02885 
##  id.1     instance     0.008665 0.09309 
##  id.2     day          0.865456 0.93030 
##  id.3     event        1.429830 1.19576 
##  Residual             10.948115 3.30879 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)  16.34097    0.91781  31.17732  17.804  < 2e-16 ***
## instance     -0.41471    0.06939 250.25579  -5.976 7.81e-09 ***
## day           1.98172    0.39305 284.93247   5.042 8.21e-07 ***
## event        -0.33323    0.24419  29.96636  -1.365    0.183    
## group        -0.62788    1.22993  28.07380  -0.511    0.614    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event 
## instance  0.136                     
## day      -0.193 -0.833              
## event    -0.123  0.011 -0.010       
## group    -0.702  0.006 -0.010 -0.010
anova(l3random_neg, l4model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_neg: negative_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
## l4model_neg: negative_score ~ 1 + instance + day + event + group + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance Chisq Df Pr(>Chisq)
## l3random_neg    9 6705.7 6751.8 -3343.9   6687.7                    
## l4model_neg    10 6707.4 6758.7 -3343.7   6687.4 0.276  1     0.5993

Explanatory Variables

Level-2 (Day)

ex1model_pos <- lmer(positive_score ~ 1 + instance + day + event + 
                       weekend + (1|id) + (0 + instance|id) + (0 + event|id), 
                     data= mlm)

summary(ex1model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + weekend + (1 |  
##     id) + (0 + instance | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 7792.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.9861 -0.5840 -0.0089  0.5645  4.0499 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 40.05787 6.3291  
##  id.1     instance     0.06107 0.2471  
##  id.2     event        7.77242 2.7879  
##  Residual             26.11287 5.1101  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   18.94303    1.27796   37.87909  14.823  < 2e-16 ***
## instance      -0.66948    0.11298  413.86975  -5.926 6.56e-09 ***
## day            3.33842    0.55108 1138.76146   6.058 1.87e-09 ***
## event          1.20518    0.52889   31.15316   2.279   0.0297 *  
## weekend       -0.01779    0.38592  319.22112  -0.046   0.9633    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event 
## instance  0.147                     
## day      -0.254 -0.869              
## event    -0.108  0.007 -0.006       
## weekend  -0.162 -0.010  0.073  0.007
anova(l3random_pos, ex1model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + event | id)
## ex1model_pos: positive_score ~ 1 + instance + day + event + weekend + (1 | id) + (0 + instance | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance Chisq Df Pr(>Chisq)
## l3random_pos    8 7807.6 7848.6 -3895.8   7791.6                    
## ex1model_pos    9 7809.6 7855.7 -3895.8   7791.6 0.003  1     0.9564
ex1model_neg <- lmer(negative_score ~ 1 + instance + day + event +
                       weekend +(1|id) + (0 + instance|id) + (0 + day|id) + 
                       (0 + event|id), data= mlm)

summary(ex1model_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + weekend + (1 |  
##     id) + (0 + instance | id) + (0 + day | id) + (0 + event |      id)
##    Data: mlm
## 
## REML criterion at convergence: 6693.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7687 -0.5528 -0.1157  0.3327  5.8376 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  8.928575 2.9881  
##  id.1     instance     0.008537 0.0924  
##  id.2     day          0.857085 0.9258  
##  id.3     event        1.431570 1.1965  
##  Residual             10.958011 3.3103  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##              Estimate Std. Error        df t value Pr(>|t|)    
## (Intercept)  15.97360    0.66272  41.27582  24.103  < 2e-16 ***
## instance     -0.41469    0.06939 251.73339  -5.976 7.77e-09 ***
## day           1.98658    0.39361 291.60959   5.047 7.91e-07 ***
## event        -0.33327    0.24432  30.04709  -1.364    0.183    
## weekend       0.06999    0.26186 336.80927   0.267    0.789    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event 
## instance  0.195                     
## day      -0.291 -0.832              
## event    -0.182  0.012 -0.009       
## weekend  -0.209 -0.006  0.063  0.010
anova(l3random_neg, ex1model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_neg: negative_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
## ex1model_neg: negative_score ~ 1 + instance + day + event + weekend + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)
## l3random_neg    9 6705.7 6751.8 -3343.9   6687.7                     
## ex1model_neg   10 6707.6 6758.9 -3343.8   6687.6 0.0785  1     0.7794

Level-3 (Event)

ex2model_pos <- lmer(positive_score ~ 1 + instance + day + event + 
                      pos_sd + pos_ac + pos_mssd + pos_pac +
                       (1|id) + (0 + instance|id) + (0 + event|id), 
                     data= mlm)

summary(ex2model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + pos_sd + pos_ac +  
##     pos_mssd + pos_pac + (1 | id) + (0 + instance | id) + (0 +      event | id)
##    Data: mlm
## 
## REML criterion at convergence: 7789.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.9518 -0.5914 -0.0123  0.5717  4.0445 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 38.84974 6.2330  
##  id.1     instance     0.05864 0.2422  
##  id.2     event        8.17175 2.8586  
##  Residual             26.08189 5.1070  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)  1.775e+01  1.710e+00  9.730e+01  10.380  < 2e-16 ***
## instance    -6.687e-01  1.126e-01  4.317e+02  -5.939 5.90e-09 ***
## day          3.321e+00  5.494e-01  1.129e+03   6.044 2.04e-09 ***
## event        5.413e-01  6.562e-01  5.996e+01   0.825   0.4127    
## pos_sd       6.196e-02  3.552e-02  1.192e+03   1.744   0.0814 .  
## pos_ac      -1.200e+01  7.327e+00  5.411e+02  -1.638   0.1021    
## pos_mssd    -3.504e-03  9.964e-03  9.860e+02  -0.352   0.7252    
## pos_pac     -1.674e+00  4.771e+00  1.041e+03  -0.351   0.7257    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event  pos_sd pos_ac ps_mss
## instance  0.112                                          
## day      -0.176 -0.873                                   
## event    -0.087  0.001  0.003                            
## pos_sd    0.062  0.010 -0.016 -0.564                     
## pos_ac    0.480  0.002  0.012  0.175 -0.272              
## pos_mssd -0.215  0.003 -0.002  0.232 -0.338  0.037       
## pos_pac  -0.280 -0.003  0.004  0.116 -0.230  0.298  0.099
anova(l3random_pos, ex2model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + event | id)
## ex2model_pos: positive_score ~ 1 + instance + day + event + pos_sd + pos_ac + pos_mssd + pos_pac + (1 | id) + (0 + instance | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)
## l3random_pos    8 7807.6 7848.6 -3895.8   7791.6                     
## ex2model_pos   12 7811.0 7872.4 -3893.5   7787.0 4.6496  4     0.3252
ex2model_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                      neg_sd + neg_ac + neg_mssd + neg_pac +
                       (1|id) + (0 + instance|id) + (0 + day|id) + (0 + event|id), 
                     data= mlm)

summary(ex2model_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + (1 | id) + (0 + instance | id) + (0 +  
##     day | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 6657.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7773 -0.5208 -0.1232  0.3182  5.7669 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  7.497391 2.73814 
##  id.1     instance     0.005555 0.07453 
##  id.2     day          0.833031 0.91271 
##  id.3     event        2.053386 1.43296 
##  Residual             10.604288 3.25642 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   12.22431    0.95836  140.02603  12.755  < 2e-16 ***
## instance      -0.41125    0.06763  277.08814  -6.080 3.96e-09 ***
## day            1.94118    0.38673  319.67575   5.019 8.62e-07 ***
## event         -1.12469    0.34894   58.44591  -3.223 0.002075 ** 
## neg_sd         0.09827    0.02944 1177.66678   3.337 0.000872 ***
## neg_ac       -16.96179    4.34401  471.07988  -3.905 0.000108 ***
## neg_mssd       0.03309    0.01136  937.10271   2.913 0.003663 ** 
## neg_pac        7.03887    2.38696  662.00438   2.949 0.003302 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event  neg_sd neg_ac ng_mss
## instance  0.131                                          
## day      -0.175 -0.841                                   
## event    -0.001  0.000  0.004                            
## neg_sd   -0.045  0.017 -0.019 -0.581                     
## neg_ac    0.650 -0.008  0.024  0.178 -0.227              
## neg_mssd -0.103  0.008  0.001  0.031  0.048  0.129       
## neg_pac  -0.405 -0.008  0.004 -0.029 -0.128 -0.050  0.036
anova(l3random_neg, ex2model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_neg: negative_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
## ex2model_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_pac + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)    
## l3random_neg    9 6705.7 6751.8 -3343.9   6687.7                         
## ex2model_neg   13 6674.9 6741.5 -3324.5   6648.9 38.811  4  7.621e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Level-4 (Participant)

ex3model_pos <- lmer(positive_score ~ 1 + instance + day + event +
                       ryff + pss + bai + bdi + mi +
                      (1|id) + (0 + instance|id) + (0 + event|id), 
                     data= mlm)

summary(ex3model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + pss + bai +  
##     bdi + mi + (1 | id) + (0 + instance | id) + (0 + event |      id)
##    Data: mlm
## 
## REML criterion at convergence: 7791.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -4.0209 -0.5799 -0.0137  0.5552  4.0907 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 34.9709  5.9136  
##  id.1     instance     0.0588  0.2425  
##  id.2     event        8.0657  2.8400  
##  Residual             26.0818  5.1070  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   -6.96581   15.33743   25.05722  -0.454   0.6536    
## instance      -0.66824    0.11260  428.16148  -5.934 6.09e-09 ***
## day            3.33047    0.54928 1138.78944   6.063 1.81e-09 ***
## event          1.20171    0.53762   31.29723   2.235   0.0327 *  
## ryff           0.31233    0.13075   24.89999   2.389   0.0248 *  
## pss           -0.12269    0.34385   24.84531  -0.357   0.7242    
## bai           -0.03522    0.22376   25.54191  -0.157   0.8762    
## bdi            0.38251    0.23530   26.23511   1.626   0.1160    
## mi             0.01790    0.19746   24.65057   0.091   0.9285    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event  ryff   pss    bai    bdi   
## instance  0.011                                                 
## day      -0.016 -0.873                                          
## event    -0.003  0.007 -0.007                                   
## ryff     -0.743  0.006 -0.009 -0.006                            
## pss      -0.587 -0.003  0.003  0.011  0.183                     
## bai       0.226  0.006 -0.006 -0.009  0.046 -0.623              
## bdi      -0.646 -0.004 -0.001 -0.016  0.545  0.214 -0.472       
## mi       -0.333 -0.004  0.003 -0.004 -0.244  0.128 -0.031  0.081
anova(l3random_pos, ex3model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + event | id)
## ex3model_pos: positive_score ~ 1 + instance + day + event + ryff + pss + bai + bdi + mi + (1 | id) + (0 + instance | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)  
## l3random_pos    8 7807.6 7848.6 -3895.8   7791.6                       
## ex3model_pos   13 7808.1 7874.7 -3891.0   7782.1 9.5108  5    0.09034 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ex4model_pos <- lmer(positive_score ~ 1 + instance + day + event +
                       ryff + 
                      (1|id) + (0 + instance|id) + (0 + event|id), 
                     data= mlm)
## Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
## Model failed to converge with max|grad| = 0.00235948 (tol = 0.002, component 1)
summary(ex4model_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + (1 | id) +  
##     (0 + instance | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 7790
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -4.0043 -0.5816 -0.0085  0.5620  4.0739 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 33.39014 5.7784  
##  id.1     instance     0.06147 0.2479  
##  id.2     event        7.85142 2.8020  
##  Residual             26.08799 5.1076  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    6.34057    5.29808   28.59371   1.197   0.2412    
## instance      -0.66744    0.11299  413.33960  -5.907 7.27e-09 ***
## day            3.33326    0.54933 1138.44935   6.068 1.76e-09 ***
## event          1.22641    0.53114   31.44420   2.309   0.0277 *  
## ryff           0.18307    0.07519   28.04986   2.435   0.0215 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event 
## instance  0.027                     
## day      -0.051 -0.870              
## event    -0.038  0.007 -0.007       
## ryff     -0.975  0.008 -0.007  0.013
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00235948 (tol = 0.002, component 1)
anova(l3random_pos, ex4model_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## l3random_pos: positive_score ~ 1 + instance + day + event + (1 | id) + (0 + instance | id) + (0 + event | id)
## ex4model_pos: positive_score ~ 1 + instance + day + event + ryff + (1 | id) + (0 + instance | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)  
## l3random_pos    8 7807.6 7848.6 -3895.8   7791.6                       
## ex4model_pos    9 7803.9 7850.0 -3892.9   7785.9 5.7206  1    0.01677 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ex3model_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                      neg_sd + neg_ac + neg_mssd + #neg_pac +
                       ryff + pss + bai + bdi + mi +
                       (1|id) + (0 + instance|id) + (0 + day|id) + (0 + event|id), 
                     data= mlm)

summary(ex3model_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + ryff + pss + bai + bdi + mi + (1 | id) + (0 +  
##     instance | id) + (0 + day | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 6677.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.7088 -0.5414 -0.1273  0.3637  5.7751 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  6.623008 2.57352 
##  id.1     instance     0.005299 0.07279 
##  id.2     day          0.685802 0.82813 
##  id.3     event        1.762343 1.32753 
##  Residual             10.741627 3.27744 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)  1.578e+01  7.387e+00  2.475e+01   2.136 0.042766 *  
## instance    -4.093e-01  6.798e-02  3.115e+02  -6.021 4.87e-09 ***
## day          1.928e+00  3.827e-01  3.583e+02   5.037 7.50e-07 ***
## event       -1.058e+00  3.362e-01  6.942e+01  -3.148 0.002420 ** 
## neg_sd       1.024e-01  2.948e-02  1.173e+03   3.473 0.000533 ***
## neg_ac      -1.424e+01  4.344e+00  4.680e+02  -3.279 0.001121 ** 
## neg_mssd     2.863e-02  1.151e-02  9.503e+02   2.488 0.013022 *  
## ryff        -3.139e-02  6.311e-02  2.480e+01  -0.497 0.623256    
## pss         -2.186e-01  1.652e-01  2.426e+01  -1.323 0.198127    
## bai          2.131e-01  1.083e-01  2.542e+01   1.968 0.060069 .  
## bdi         -3.275e-03  1.155e-01  2.754e+01  -0.028 0.977576    
## mi           6.666e-02  9.480e-02  2.415e+01   0.703 0.488647    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event  neg_sd neg_ac ng_mss ryff   pss    bai   
## instance  0.010                                                               
## day      -0.011 -0.856                                                        
## event    -0.014 -0.001  0.004                                                 
## neg_sd    0.028  0.015 -0.017 -0.613                                          
## neg_ac    0.038 -0.008  0.022  0.182 -0.243                                   
## neg_mssd  0.090  0.007  0.003  0.028  0.067  0.102                            
## ryff     -0.744  0.012 -0.015  0.012 -0.042  0.033 -0.103                     
## pss      -0.577 -0.003  0.003 -0.013  0.041 -0.082  0.042  0.171              
## bai       0.223  0.006 -0.005  0.015 -0.042  0.073 -0.059  0.055 -0.629       
## bdi      -0.647  0.003 -0.012  0.015 -0.063  0.089 -0.123  0.554  0.197 -0.446
## mi       -0.331 -0.001  0.000  0.011 -0.027  0.064 -0.044 -0.234  0.116 -0.024
##          bdi   
## instance       
## day            
## event          
## neg_sd         
## neg_ac         
## neg_mssd       
## ryff           
## pss            
## bai            
## bdi            
## mi        0.095
anova(ex2model_neg, ex3model_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## ex2model_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_pac + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
## ex3model_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + ryff + pss + bai + bdi + mi + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##              npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)
## ex2model_neg   13 6674.9 6741.5 -3324.5   6648.9                     
## ex3model_neg   17 6682.3 6769.4 -3324.1   6648.3 0.6302  4     0.9597

Cross-Level Interactions

Treatment Group

crosslevel1_pos <- lmer(positive_score ~ 1 + instance + day + event + ryff +
                          group + group:instance + (1|id), data = mlm)

summary(crosslevel1_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + group +  
##     group:instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 7971.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3442 -0.6189 -0.0472  0.5166  4.3267 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 28.13    5.304   
##  Residual             33.12    5.755   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                  Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)      17.86179    4.76478   30.10820   3.749 0.000755 ***
## instance         -0.70755    0.12330 1204.73081  -5.739 1.21e-08 ***
## day               3.27851    0.61643 1204.77173   5.319 1.25e-07 ***
## event             1.37024    0.20409 1204.97446   6.714 2.91e-11 ***
## ryff              0.01772    0.06463   28.80623   0.274 0.785932    
## group            -0.64520    2.00677   34.92727  -0.322 0.749739    
## instance:group    0.08330    0.07615 1204.53223   1.094 0.274221    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  ryff   group 
## instance     0.025                                   
## day         -0.074 -0.897                            
## event       -0.087  0.014 -0.014                     
## ryff        -0.945 -0.002  0.002  0.000              
## group       -0.286  0.104 -0.008  0.002  0.069       
## instanc:grp  0.067 -0.324  0.009  0.013 -0.003 -0.301
crosslevel1_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          group + group:instance + (1|id), data = mlm)

summary(crosslevel1_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + group + group:instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6752.6
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9132 -0.5496 -0.1388  0.3736  6.2479 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.28    2.878   
##  Residual             12.43    3.525   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                  Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)      14.33722    1.05730   98.63297  13.560  < 2e-16 ***
## instance         -0.43059    0.07555 1201.03107  -5.699 1.51e-08 ***
## day               1.88615    0.37778 1201.02009   4.993 6.83e-07 ***
## event            -0.83364    0.24788 1213.62927  -3.363 0.000795 ***
## neg_sd            0.09240    0.03090 1217.73177   2.990 0.002847 ** 
## neg_ac           -9.18387    3.53414  880.34531  -2.599 0.009517 ** 
## neg_mssd          0.03318    0.01074 1229.89580   3.089 0.002052 ** 
## neg_pac           1.90937    2.01447 1229.94735   0.948 0.343404    
## group            -0.78503    1.10765   38.51991  -0.709 0.482754    
## instance:group    0.04803    0.04665 1201.33440   1.030 0.303450    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc group 
## instance     0.053                                                        
## day         -0.185 -0.897                                                 
## event       -0.110 -0.011  0.011                                          
## neg_sd       0.017  0.021 -0.021 -0.856                                   
## neg_ac       0.442 -0.027  0.031  0.167 -0.181                            
## neg_mssd    -0.109  0.004 -0.003  0.047  0.021  0.122                     
## neg_pac     -0.329  0.002  0.001  0.011 -0.123 -0.016  0.036              
## group       -0.534  0.113 -0.007  0.014 -0.018  0.070  0.022  0.042       
## instanc:grp  0.181 -0.324  0.009  0.012 -0.004  0.010  0.014 -0.010 -0.333
crosslevel2_pos <- lmer(positive_score ~ 1 + instance + day + event + ryff +
                          group + group:day + (1|id), data = mlm)

summary(crosslevel2_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + group +  
##     group:day + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 7969
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3320 -0.6206 -0.0507  0.5165  4.3203 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 28.12    5.303   
##  Residual             33.15    5.757   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   17.69957    4.77148   30.29871   3.709 0.000834 ***
## instance      -0.66419    0.11668 1205.03806  -5.692 1.57e-08 ***
## day            3.18260    0.64860 1204.80586   4.907 1.05e-06 ***
## event          1.36842    0.20417 1204.97351   6.702 3.14e-11 ***
## ryff           0.01787    0.06462   28.80890   0.277 0.784102    
## group         -0.34499    2.07517   39.94493  -0.166 0.868801    
## day:group      0.18056    0.40261 1204.29072   0.448 0.653891    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##           (Intr) instnc day    event  ryff   group 
## instance   0.049                                   
## day       -0.097 -0.896                            
## event     -0.086  0.019 -0.017                     
## ryff      -0.944 -0.003  0.003  0.000              
## group     -0.290  0.009  0.115  0.002  0.066       
## day:group  0.087 -0.007 -0.310  0.011 -0.002 -0.387
crosslevel2_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          group + group:day + (1|id), data = mlm)

summary(crosslevel2_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + group + group:day + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6748.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9048 -0.5426 -0.1360  0.3758  6.2279 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.279   2.877   
##  Residual             12.421   3.524   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   14.45488    1.07123  103.66376  13.494  < 2e-16 ***
## instance      -0.40603    0.07146 1201.19247  -5.682 1.67e-08 ***
## day            1.73275    0.39719 1201.03811   4.363 1.40e-05 ***
## event         -0.83417    0.24783 1213.61463  -3.366 0.000787 ***
## neg_sd         0.09249    0.03090 1217.72147   2.993 0.002814 ** 
## neg_ac        -9.12873    3.53417  880.08375  -2.583 0.009955 ** 
## neg_mssd       0.03323    0.01074 1229.89447   3.094 0.002021 ** 
## neg_pac        1.90768    2.01411 1229.94818   0.947 0.343743    
## group         -1.00299    1.15353   45.28943  -0.869 0.389160    
## day:group      0.30080    0.24653 1201.19944   1.220 0.222649    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##           (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc group 
## instance   0.115                                                        
## day       -0.249 -0.896                                                 
## event     -0.109 -0.007  0.008                                          
## neg_sd     0.017  0.020 -0.020 -0.856                                   
## neg_ac     0.440 -0.025  0.023  0.167 -0.181                            
## neg_mssd  -0.107  0.009 -0.007  0.047  0.021  0.122                     
## neg_pac   -0.325 -0.002  0.004  0.011 -0.123 -0.016  0.036              
## group     -0.551  0.009  0.128  0.014 -0.018  0.061  0.019  0.041       
## day:group  0.241 -0.008 -0.309  0.008 -0.001  0.021  0.015 -0.009 -0.425
crosslevel3_pos <- lmer(positive_score ~ 1 + instance + day + event + ryff +
                          group + group:event + (1|id), data = mlm)

summary(crosslevel3_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + group +  
##     group:event + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 7968.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3536 -0.6190 -0.0410  0.5227  4.2892 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 28.12    5.303   
##  Residual             33.13    5.756   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   17.92750    4.77337   30.34532   3.756 0.000734 ***
## instance      -0.66290    0.11665 1205.05041  -5.683 1.66e-08 ***
## day            3.27036    0.61649 1204.78595   5.305 1.34e-07 ***
## event          1.17171    0.28969 1204.31173   4.045 5.57e-05 ***
## ryff           0.01766    0.06462   28.81077   0.273 0.786565    
## group         -0.76251    2.08063   40.35943  -0.366 0.715920    
## event:group    0.38855    0.40814 1204.91144   0.952 0.341290    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  ryff   group 
## instance     0.050                                   
## day         -0.074 -0.945                            
## event       -0.126  0.008 -0.007                     
## ryff        -0.944 -0.003  0.002  0.003              
## group       -0.292  0.003 -0.004  0.283  0.067       
## event:group  0.091  0.008 -0.004 -0.710 -0.004 -0.393
crosslevel3_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          group + group:event + (1|id), data = mlm)

summary(crosslevel3_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + group + group:event + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6750.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9553 -0.5384 -0.1460  0.3751  6.2792 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.285   2.878   
##  Residual             12.435   3.526   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   14.22749    1.07464  104.33635  13.239  < 2e-16 ***
## instance      -0.40518    0.07150 1201.21063  -5.667 1.82e-08 ***
## day            1.88240    0.37791 1201.02357   4.981 7.25e-07 ***
## event         -0.87719    0.27754 1209.18120  -3.161  0.00161 ** 
## neg_sd         0.09251    0.03092 1217.73546   2.992  0.00283 ** 
## neg_ac        -9.12545    3.54702  877.20701  -2.573  0.01025 *  
## neg_mssd       0.03347    0.01083 1229.89869   3.090  0.00204 ** 
## neg_pac        1.94066    2.01543 1229.95606   0.963  0.33578    
## group         -0.56596    1.15635   45.66281  -0.489  0.62687    
## event:group    0.08194    0.25258 1203.60702   0.324  0.74568    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc group 
## instance     0.118                                                        
## day         -0.184 -0.945                                                 
## event       -0.212 -0.010  0.011                                          
## neg_sd       0.016  0.020 -0.021 -0.763                                   
## neg_ac       0.453 -0.024  0.030  0.112 -0.180                            
## neg_mssd    -0.078  0.010 -0.003 -0.014  0.021  0.130                     
## neg_pac     -0.318 -0.002  0.001  0.003 -0.123 -0.015  0.037              
## group       -0.554  0.002 -0.003  0.208 -0.017  0.035 -0.028  0.030       
## event:group  0.251  0.008 -0.002 -0.449 -0.004  0.081  0.124  0.016 -0.429

Significant Predictors

crosslevel4_pos <- lmer(positive_score ~ 1 + instance + day + event + 
                          ryff + ryff:instance + (1|id), data = mlm)

summary(crosslevel4_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + ryff:instance +  
##     (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 7981.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3000 -0.6242 -0.0419  0.5204  4.3543 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 27.20    5.215   
##  Residual             33.13    5.756   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                 Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    1.634e+01  4.700e+00  3.721e+01   3.476  0.00131 ** 
## instance      -5.087e-01  2.154e-01  1.205e+03  -2.362  0.01834 *  
## day            3.260e+00  6.167e-01  1.205e+03   5.286 1.48e-07 ***
## event          1.367e+00  2.041e-01  1.205e+03   6.696 3.27e-11 ***
## ryff           3.516e-02  6.654e-02  3.609e+01   0.528  0.60043    
## instance:ryff -2.200e-03  2.566e-03  1.204e+03  -0.857  0.39139    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  ryff  
## instance    -0.219                            
## day         -0.069 -0.533                     
## event       -0.086  0.007 -0.014              
## ryff        -0.971  0.252 -0.005 -0.001       
## instnc:ryff  0.294 -0.841  0.025  0.004 -0.302
crosslevel5_pos <- lmer(positive_score ~ 1 + instance + day + event + 
                          ryff + ryff:day + (1|id), data = mlm)

summary(crosslevel5_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + ryff:day +  
##     (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 7979.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3157 -0.6185 -0.0478  0.5103  4.3248 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 27.15    5.210   
##  Residual             33.15    5.758   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)  1.757e+01  4.849e+00  4.227e+01   3.624 0.000772 ***
## instance    -6.641e-01  1.167e-01  1.205e+03  -5.690 1.60e-08 ***
## day          3.248e+00  1.107e+00  1.204e+03   2.934 0.003405 ** 
## event        1.368e+00  2.042e-01  1.205e+03   6.698 3.24e-11 ***
## ryff         1.718e-02  6.875e-02  4.128e+01   0.250 0.803938    
## day:ryff     3.776e-04  1.349e-02  1.204e+03   0.028 0.977672    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr) instnc day    event  ryff  
## instance  0.041                            
## day      -0.356 -0.506                     
## event    -0.081  0.019 -0.015              
## ryff     -0.973  0.006  0.323 -0.003       
## day:ryff  0.378 -0.025 -0.830  0.008 -0.388
crosslevel6_pos <- lmer(positive_score ~ 1 + instance + day + event + 
                          ryff + ryff:event + (1|id), data = mlm)

summary(crosslevel6_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + ryff:event +  
##     (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 7934.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.1740 -0.6409 -0.0446  0.5354  4.0126 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 27.24    5.219   
##  Residual             31.95    5.652   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    4.89394    4.86611   42.70454   1.006  0.32022    
## instance      -0.65463    0.11455 1205.14231  -5.715 1.38e-08 ***
## day            3.21549    0.60543 1204.86569   5.311 1.30e-07 ***
## event          7.64697    0.95286 1204.86531   8.025 2.39e-15 ***
## ryff           0.20032    0.06896   41.61795   2.905  0.00586 ** 
## event:ryff    -0.09039    0.01341 1204.66271  -6.741 2.43e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##            (Intr) instnc day    event  ryff  
## instance    0.045                            
## day        -0.068 -0.945                     
## event      -0.394  0.016 -0.017              
## ryff       -0.973  0.002 -0.003  0.384       
## event:ryff  0.385 -0.012  0.014 -0.978 -0.392
intmodel1_pos <- lmer(positive_score ~ 1 + instance + day + event +
                       ryff + ryff:event + (1|id) + (0 + instance|id) + (0 + event|id), 
                     data= mlm)

summary(intmodel1_pos)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: positive_score ~ 1 + instance + day + event + ryff + ryff:event +  
##     (1 | id) + (0 + instance | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 7788.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.9956 -0.5836 -0.0055  0.5517  4.1257 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept) 31.87258 5.6456  
##  id.1     instance     0.06167 0.2483  
##  id.2     event        6.38992 2.5278  
##  Residual             26.12131 5.1109  
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    4.80321    5.23257   29.44857   0.918  0.36611    
## instance      -0.66535    0.11307  412.48288  -5.885 8.27e-09 ***
## day            3.32080    0.54965 1136.03257   6.042 2.06e-09 ***
## event          7.07678    2.30569   29.93978   3.069  0.00453 ** 
## ryff           0.20480    0.07425   28.89020   2.758  0.00998 ** 
## event:ryff    -0.08491    0.03275   29.55394  -2.592  0.01467 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##            (Intr) instnc day    event  ryff  
## instance    0.027                            
## day        -0.051 -0.870                     
## event      -0.135  0.007 -0.008              
## ryff       -0.975  0.009 -0.008  0.130       
## event:ryff  0.130 -0.005  0.007 -0.978 -0.131
anova(ex4model_pos,intmodel1_pos)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## ex4model_pos: positive_score ~ 1 + instance + day + event + ryff + (1 | id) + (0 + instance | id) + (0 + event | id)
## intmodel1_pos: positive_score ~ 1 + instance + day + event + ryff + ryff:event + (1 | id) + (0 + instance | id) + (0 + event | id)
##               npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)  
## ex4model_pos     9 7803.9 7850.0 -3892.9   7785.9                       
## intmodel1_pos   10 7799.4 7850.6 -3889.7   7779.4 6.5174  1    0.01068 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
crosslevel4_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_sd:instance + (1|id), data = mlm)
## Warning: Some predictor variables are on very different scales: consider
## rescaling

## Warning: Some predictor variables are on very different scales: consider
## rescaling
summary(crosslevel4_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_sd:instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6757.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9945 -0.5305 -0.1410  0.3768  6.2211 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.065   2.840   
##  Residual             12.404   3.522   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                   Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)      1.337e+01  9.475e-01  2.431e+02  14.109  < 2e-16 ***
## instance        -3.336e-01  8.243e-02  1.201e+03  -4.047 5.52e-05 ***
## day              1.886e+00  3.774e-01  1.201e+03   4.998 6.66e-07 ***
## event           -8.353e-01  2.476e-01  1.214e+03  -3.374 0.000765 ***
## neg_sd           1.416e-01  4.187e-02  1.211e+03   3.381 0.000744 ***
## neg_ac          -8.968e+00  3.516e+00  8.572e+02  -2.551 0.010922 *  
## neg_mssd         3.326e-02  1.072e-02  1.231e+03   3.101 0.001971 ** 
## neg_pac          1.961e+00  2.010e+00  1.231e+03   0.976 0.329428    
## instance:neg_sd -6.130e-03  3.526e-03  1.201e+03  -1.738 0.082388 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance    -0.055                                                 
## day         -0.212 -0.816                                          
## event       -0.115 -0.006  0.011                                   
## neg_sd      -0.227  0.351 -0.011 -0.630                            
## neg_ac       0.528 -0.011  0.031  0.165 -0.117                     
## neg_mssd    -0.111  0.011 -0.003  0.047  0.020  0.120              
## neg_pac     -0.342 -0.001  0.001  0.011 -0.090 -0.019  0.035       
## instnc:ng_s  0.345 -0.500 -0.007 -0.001 -0.676 -0.022 -0.006 -0.001
## fit warnings:
## Some predictor variables are on very different scales: consider rescaling
crosslevel5_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_sd:day + (1|id), data = mlm)

summary(crosslevel5_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_sd:day + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6752.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.0148 -0.5347 -0.1481  0.3765  6.1561 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.067   2.840   
##  Residual             12.382   3.519   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   12.94489    0.99041  282.63139  13.070  < 2e-16 ***
## instance      -0.40491    0.07134 1201.27954  -5.676 1.73e-08 ***
## day            2.38023    0.43632 1200.77361   5.455 5.93e-08 ***
## event         -0.83637    0.24737 1214.39172  -3.381 0.000745 ***
## neg_sd         0.17735    0.04849 1208.58268   3.657 0.000266 ***
## neg_ac        -9.01305    3.51253  858.01145  -2.566 0.010458 *  
## neg_mssd       0.03327    0.01071 1230.96723   3.105 0.001944 ** 
## neg_pac        1.95656    2.00858 1230.67756   0.974 0.330198    
## day:neg_sd    -0.04227    0.01862 1200.63587  -2.270 0.023368 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##            (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance    0.128                                                 
## day        -0.395 -0.816                                          
## event      -0.109 -0.007  0.009                                   
## neg_sd     -0.335  0.014  0.376 -0.545                            
## neg_ac      0.507 -0.026  0.032  0.165 -0.105                     
## neg_mssd   -0.106  0.009  0.001  0.047  0.018  0.120              
## neg_pac    -0.326 -0.002  0.001  0.011 -0.078 -0.019  0.035       
## day:neg_sd  0.441 -0.002 -0.503  0.002 -0.772 -0.011 -0.006  0.001
crosslevel6_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_ac:instance + (1|id), data = mlm)

summary(crosslevel6_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_ac:instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6748.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9610 -0.5424 -0.1498  0.3807  6.2717 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.031   2.834   
##  Residual             12.405   3.522   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                   Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)       14.82891    1.02563  306.58330  14.458  < 2e-16 ***
## instance          -0.52799    0.10036 1201.74484  -5.261 1.69e-07 ***
## day                1.89774    0.37756 1201.09435   5.026 5.76e-07 ***
## event             -0.82884    0.24762 1214.58072  -3.347 0.000841 ***
## neg_sd             0.09246    0.03087 1218.68948   2.995 0.002797 ** 
## neg_ac            -2.88513    5.00744 1111.88974  -0.576 0.564618    
## neg_mssd           0.03318    0.01072 1230.97688   3.094 0.002016 ** 
## neg_pac            1.98685    2.01036 1230.63523   0.988 0.323195    
## instance:neg_ac   -0.81874    0.47008 1203.51713  -1.742 0.081819 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance    -0.262                                                 
## day         -0.182 -0.689                                          
## event       -0.098 -0.016  0.012                                   
## neg_sd       0.008  0.014 -0.022 -0.855                            
## neg_ac       0.703 -0.513  0.039  0.127 -0.125                     
## neg_mssd    -0.099  0.005 -0.003  0.047  0.022  0.086              
## neg_pac     -0.311 -0.007  0.002  0.011 -0.122 -0.008  0.035       
## instnc:ng_c -0.499  0.703 -0.024 -0.016 -0.001 -0.712 -0.003 -0.008
crosslevel7_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_ac:day + (1|id), data = mlm)

summary(crosslevel7_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_ac:day + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6746.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9527 -0.5406 -0.1485  0.3829  6.2646 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.045   2.836   
##  Residual             12.427   3.525   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   14.66942    1.18909  463.27508  12.337  < 2e-16 ***
## instance      -0.40788    0.07153 1201.28882  -5.702 1.49e-08 ***
## day            1.53126    0.53438 1203.11613   2.865  0.00424 ** 
## event         -0.83061    0.24787 1214.66413  -3.351  0.00083 ***
## neg_sd         0.09199    0.03090 1218.84541   2.977  0.00297 ** 
## neg_ac        -4.00973    6.51479 1194.34348  -0.615  0.53835    
## neg_mssd       0.03323    0.01073 1230.97349   3.096  0.00200 ** 
## neg_pac        1.95950    2.01206 1230.63531   0.974  0.33031    
## day:neg_ac    -2.50378    2.69836 1205.79657  -0.928  0.35365    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##            (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance    0.080                                                 
## day        -0.588 -0.639                                          
## event      -0.077 -0.008 -0.007                                   
## neg_sd     -0.003  0.021 -0.005 -0.856                            
## neg_ac      0.789 -0.048 -0.583  0.108 -0.109                     
## neg_mssd   -0.080  0.009 -0.009  0.047  0.022  0.073              
## neg_pac    -0.272 -0.002  0.000  0.011 -0.122 -0.010  0.035       
## day:neg_ac -0.664  0.041  0.707 -0.022  0.015 -0.842 -0.010 -0.001
crosslevel8_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_mssd:instance + (1|id), data = mlm)
## Warning: Some predictor variables are on very different scales: consider
## rescaling

## Warning: Some predictor variables are on very different scales: consider
## rescaling
summary(crosslevel8_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_mssd:instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6755.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9692 -0.5442 -0.1443  0.3684  6.3200 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.023   2.832   
##  Residual             12.372   3.517   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                     Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)        1.446e+01  9.119e-01  2.074e+02  15.859  < 2e-16 ***
## instance          -4.584e-01  7.438e-02  1.201e+03  -6.163 9.71e-10 ***
## day                1.889e+00  3.770e-01  1.201e+03   5.010 6.25e-07 ***
## event             -8.341e-01  2.473e-01  1.214e+03  -3.373 0.000766 ***
## neg_sd             9.224e-02  3.083e-02  1.219e+03   2.992 0.002825 ** 
## neg_ac            -8.459e+00  3.519e+00  8.436e+02  -2.404 0.016443 *  
## neg_mssd          -4.006e-03  1.822e-02  1.224e+03  -0.220 0.826064    
## neg_pac            1.840e+00  2.008e+00  1.231e+03   0.916 0.359593    
## instance:neg_mssd  4.865e-03  1.931e-03  1.210e+03   2.519 0.011898 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance     0.069                                                 
## day         -0.216 -0.908                                          
## event       -0.118 -0.008  0.011                                   
## neg_sd       0.008  0.020 -0.022 -0.855                            
## neg_ac       0.571 -0.045  0.032  0.165 -0.179                     
## neg_mssd    -0.251  0.235 -0.007  0.025  0.015  0.012              
## neg_pac     -0.360  0.005  0.001  0.010 -0.122 -0.021  0.039       
## instnc:ng_m  0.228 -0.284  0.007  0.002 -0.002  0.072 -0.809 -0.023
## fit warnings:
## Some predictor variables are on very different scales: consider rescaling
crosslevel9_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_mssd:day + (1|id), data = mlm)

summary(crosslevel9_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_mssd:day + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6755.6
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9609 -0.5417 -0.1470  0.3692  6.3187 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.03    2.834   
##  Residual             12.40    3.522   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   1.439e+01  9.248e-01  2.176e+02  15.560  < 2e-16 ***
## instance     -4.059e-01  7.141e-02  1.201e+03  -5.685 1.64e-08 ***
## day           1.691e+00  3.926e-01  1.202e+03   4.306 1.79e-05 ***
## event        -8.363e-01  2.476e-01  1.214e+03  -3.378 0.000754 ***
## neg_sd        9.242e-02  3.087e-02  1.219e+03   2.994 0.002809 ** 
## neg_ac       -8.742e+00  3.520e+00  8.472e+02  -2.484 0.013198 *  
## neg_mssd     -2.457e-03  2.281e-02  1.218e+03  -0.108 0.914249    
## neg_pac       1.893e+00  2.011e+00  1.231e+03   0.942 0.346537    
## day:neg_mssd  1.830e-02  1.035e-02  1.207e+03   1.768 0.077363 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance     0.137                                                 
## day         -0.283 -0.907                                          
## event       -0.118 -0.007  0.011                                   
## neg_sd       0.009  0.020 -0.021 -0.855                            
## neg_ac       0.563 -0.026  0.014  0.165 -0.179                     
## neg_mssd    -0.297  0.009  0.242  0.023  0.010  0.006              
## neg_pac     -0.355 -0.002  0.006  0.011 -0.122 -0.020  0.032       
## day:ng_mssd  0.277 -0.006 -0.275 -0.002  0.000  0.057 -0.883 -0.018
crosslevel10_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_pac:instance + (1|id), data = mlm)

summary(crosslevel10_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_pac:instance + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6745.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.0067 -0.5342 -0.1447  0.3891  6.2623 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.055   2.838   
##  Residual             12.371   3.517   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                    Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)        12.76085    1.00420  299.12864  12.708  < 2e-16 ***
## instance           -0.25151    0.09396 1201.28698  -2.677 0.007537 ** 
## day                 1.87211    0.37695 1201.08325   4.966  7.8e-07 ***
## event              -0.83862    0.24726 1214.38863  -3.392 0.000717 ***
## neg_sd              0.09288    0.03083 1218.61274   3.013 0.002642 ** 
## neg_ac             -8.94599    3.51112  857.43492  -2.548 0.011011 *  
## neg_mssd            0.03339    0.01071 1230.96733   3.117 0.001867 ** 
## neg_pac             8.97535    3.44047 1217.79357   2.609 0.009198 ** 
## instance:neg_pac   -0.88592    0.35278 1200.75902  -2.511 0.012160 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance    -0.207                                                 
## day         -0.193 -0.724                                          
## event       -0.106 -0.009  0.011                                   
## neg_sd       0.005  0.019 -0.022 -0.855                            
## neg_ac       0.496 -0.008  0.031  0.165 -0.179                     
## neg_mssd    -0.107  0.013 -0.003  0.047  0.022  0.120              
## neg_pac     -0.566  0.528 -0.008  0.002 -0.066  0.003  0.028       
## instnc:ng_p  0.466 -0.651  0.010  0.005 -0.006 -0.018 -0.009 -0.812
crosslevel11_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                          neg_sd + neg_ac + neg_mssd + neg_pac +
                          neg_pac:day + (1|id), data = mlm)

summary(crosslevel11_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_pac + neg_pac:day + (1 | id)
##    Data: mlm
## 
## REML criterion at convergence: 6742.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.0103 -0.5337 -0.1496  0.3860  6.2470 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  id       (Intercept)  8.053   2.838   
##  Residual             12.374   3.518   
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##               Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   12.37688    1.09237  390.58501  11.330  < 2e-16 ***
## instance      -0.40156    0.07133 1201.29129  -5.629 2.25e-08 ***
## day            2.65656    0.49167 1200.81896   5.403 7.89e-08 ***
## event         -0.83809    0.24728 1214.39514  -3.389 0.000724 ***
## neg_sd         0.09282    0.03083 1218.62041   3.011 0.002661 ** 
## neg_ac        -9.00561    3.51113  857.39198  -2.565 0.010490 *  
## neg_mssd       0.03336    0.01071 1230.96830   3.115 0.001884 ** 
## neg_pac       11.20977    4.27107 1212.26296   2.625 0.008785 ** 
## day:neg_pac   -4.64179    1.89150 1200.67696  -2.454 0.014268 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss neg_pc
## instance     0.105                                                 
## day         -0.513 -0.711                                          
## event       -0.097 -0.007  0.006                                   
## neg_sd       0.004  0.020 -0.013 -0.855                            
## neg_ac       0.458 -0.025  0.031  0.165 -0.179                     
## neg_mssd    -0.099  0.009  0.003  0.047  0.022  0.120              
## neg_pac     -0.653  0.017  0.567  0.001 -0.053  0.001  0.024       
## day:neg_pac  0.582 -0.021 -0.642  0.004 -0.006 -0.011 -0.008 -0.883
#add significant interactions
intmodel1_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                      neg_sd + neg_ac + neg_mssd + neg_sd:day + neg_mssd:day +
                        + neg_pac:instance + neg_pac:day + (1|id) + 
                        (0 + instance|id) + (0 + day|id) + (0 + event|id), 
                     data= mlm)

summary(intmodel1_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_sd:day + neg_mssd:day + +neg_pac:instance +  
##     neg_pac:day + (1 | id) + (0 + instance | id) + (0 + day |  
##     id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 6665.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.8702 -0.5089 -0.1307  0.3504  5.6821 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  7.045189 2.65428 
##  id.1     instance     0.002467 0.04967 
##  id.2     day          0.847794 0.92076 
##  id.3     event        1.850278 1.36025 
##  Residual             10.632598 3.26077 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                    Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)       1.290e+01  1.007e+00  1.766e+02  12.812  < 2e-16 ***
## instance         -1.021e-01  1.727e-01  5.514e+02  -0.591 0.554700    
## day               6.660e-01  8.392e-01  9.606e+02   0.794 0.427648    
## event            -1.108e+00  3.396e-01  6.425e+01  -3.263 0.001769 ** 
## neg_sd            1.834e-01  4.594e-02  1.175e+03   3.992 6.96e-05 ***
## neg_ac           -1.587e+01  4.312e+00  4.472e+02  -3.680 0.000262 ***
## neg_mssd          6.991e-04  2.365e-02  8.476e+02   0.030 0.976421    
## day:neg_sd       -3.982e-02  1.772e-02  1.152e+03  -2.247 0.024805 *  
## day:neg_mssd      1.699e-02  1.098e-02  7.610e+02   1.548 0.122002    
## instance:neg_pac -1.829e+00  9.375e-01  6.550e+02  -1.951 0.051531 .  
## day:neg_pac       9.274e+00  4.110e+00  8.649e+02   2.257 0.024282 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss dy:ng_s dy:ng_m
## instance     0.079                                                          
## day         -0.251 -0.902                                                   
## event       -0.013 -0.012  0.015                                            
## neg_sd      -0.395 -0.064  0.254 -0.387                                     
## neg_ac       0.614  0.005  0.003  0.180 -0.152                              
## neg_mssd    -0.308  0.034  0.113  0.015  0.072  0.003                       
## day:neg_sd   0.440  0.074 -0.306  0.005 -0.768  0.011 -0.082                
## day:ng_mssd  0.306 -0.026 -0.142  0.000 -0.066  0.068 -0.878  0.094         
## instnc:ng_p -0.036 -0.922  0.816  0.013  0.075 -0.009 -0.032 -0.080   0.023 
## day:neg_pac  0.027  0.890 -0.836 -0.018 -0.071 -0.004  0.034  0.053  -0.023 
##             inst:_
## instance          
## day               
## event             
## neg_sd            
## neg_ac            
## neg_mssd          
## day:neg_sd        
## day:ng_mssd       
## instnc:ng_p       
## day:neg_pac -0.967
anova(ex2model_neg, intmodel1_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## ex2model_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_pac + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
## intmodel1_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_sd:day + neg_mssd:day + +neg_pac:instance + neg_pac:day + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##               npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)
## ex2model_neg    13 6674.9 6741.5 -3324.5   6648.9                     
## intmodel1_neg   16 6674.8 6756.7 -3321.4   6642.8 6.1576  3     0.1042
#only significant
intmodel2_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                      neg_sd + neg_ac + neg_mssd + 
                        neg_sd:day + neg_pac:instance + neg_pac:day +
                       (1|id) + (0 + instance|id) + (0 + day|id) + (0 + event|id), 
                     data= mlm)

summary(intmodel2_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_sd:day + neg_pac:instance + neg_pac:day +  
##     (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event |      id)
##    Data: mlm
## 
## REML criterion at convergence: 6660.6
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.9078 -0.5212 -0.1372  0.3475  5.6861 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  7.247271 2.69208 
##  id.1     instance     0.004655 0.06823 
##  id.2     day          0.829043 0.91052 
##  id.3     event        1.873975 1.36893 
##  Residual             10.625327 3.25965 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                    Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)        12.43608    0.96325  160.03732  12.910  < 2e-16 ***
## instance           -0.09378    0.17349  569.11782  -0.541  0.58904    
## day                 0.84287    0.83238 1010.54874   1.013  0.31150    
## event              -1.10895    0.34067   63.43081  -3.255  0.00182 ** 
## neg_sd              0.18778    0.04584 1168.36049   4.096 4.49e-05 ***
## neg_ac            -16.28637    4.31213  456.85357  -3.777  0.00018 ***
## neg_mssd            0.03287    0.01133  931.37766   2.902  0.00379 ** 
## day:neg_sd         -0.04224    0.01764 1141.98377  -2.395  0.01677 *  
## instance:neg_pac   -1.86726    0.94112  683.31014  -1.984  0.04765 *  
## day:neg_pac         9.44004    4.12154  891.60207   2.290  0.02223 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss dy:ng_s inst:_
## instance     0.091                                                         
## day         -0.219 -0.915                                                  
## event       -0.013 -0.012  0.016                                           
## neg_sd      -0.393 -0.066  0.248 -0.387                                    
## neg_ac       0.623  0.007  0.012  0.180 -0.149                             
## neg_mssd    -0.086  0.025 -0.024  0.032  0.031  0.131                      
## day:neg_sd   0.432  0.077 -0.297  0.006 -0.767  0.005  0.001               
## instnc:ng_p -0.045 -0.921  0.829  0.013  0.077 -0.011 -0.024 -0.083        
## day:neg_pac  0.036  0.890 -0.849 -0.018 -0.073 -0.002  0.028  0.056  -0.967
anova(ex2model_neg, intmodel2_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## ex2model_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_pac + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
## intmodel2_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_sd:day + neg_pac:instance + neg_pac:day + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##               npar    AIC    BIC  logLik deviance  Chisq Df Pr(>Chisq)
## ex2model_neg    13 6674.9 6741.5 -3324.5   6648.9                     
## intmodel2_neg   15 6675.2 6752.0 -3322.6   6645.2 3.7046  2     0.1569
#only mssd
intmodel3_neg <- lmer(negative_score ~ 1 + instance + day + event + 
                      neg_sd + neg_ac + neg_mssd + neg_mssd:day +
                       (1|id) + (0 + instance|id) + (0 + day|id) + (0 + event|id), 
                     data= mlm)

summary(intmodel3_neg)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac +  
##     neg_mssd + neg_mssd:day + (1 | id) + (0 + instance | id) +  
##     (0 + day | id) + (0 + event | id)
##    Data: mlm
## 
## REML criterion at convergence: 6673.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.6829 -0.5233 -0.1181  0.3557  5.7594 
## 
## Random effects:
##  Groups   Name        Variance  Std.Dev.
##  id       (Intercept)  7.291017 2.70019 
##  id.1     instance     0.002552 0.05051 
##  id.2     day          0.825339 0.90848 
##  id.3     event        1.661083 1.28883 
##  Residual             10.728271 3.27540 
## Number of obs: 1240, groups:  id, 32
## 
## Fixed effects:
##                Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)   1.393e+01  9.081e-01  1.216e+02  15.339  < 2e-16 ***
## instance     -4.112e-01  6.728e-02  2.972e+02  -6.111 3.10e-09 ***
## day           1.727e+00  4.048e-01  3.722e+02   4.265 2.54e-05 ***
## event        -1.073e+00  3.311e-01  7.229e+01  -3.242 0.001797 ** 
## neg_sd        1.087e-01  2.931e-02  1.182e+03   3.708 0.000219 ***
## neg_ac       -1.517e+01  4.289e+00  4.390e+02  -3.536 0.000449 ***
## neg_mssd     -5.485e-03  2.364e-02  8.550e+02  -0.232 0.816553    
## day:neg_mssd  1.992e-02  1.096e-02  7.662e+02   1.818 0.069417 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) instnc day    event  neg_sd neg_ac ng_mss
## instance     0.132                                          
## day         -0.260 -0.813                                   
## event       -0.020 -0.001  0.005                            
## neg_sd      -0.095  0.016 -0.021 -0.622                     
## neg_ac       0.675 -0.011  0.006  0.182 -0.230              
## neg_mssd    -0.304  0.013  0.249  0.016  0.017  0.003       
## day:ng_mssd  0.295 -0.011 -0.284 -0.001  0.010  0.068 -0.878
anova(ex2model_neg, intmodel3_neg)
## refitting model(s) with ML (instead of REML)
## Data: mlm
## Models:
## ex2model_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_pac + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
## intmodel3_neg: negative_score ~ 1 + instance + day + event + neg_sd + neg_ac + neg_mssd + neg_mssd:day + (1 | id) + (0 + instance | id) + (0 + day | id) + (0 + event | id)
##               npar    AIC    BIC  logLik deviance Chisq Df Pr(>Chisq)
## ex2model_neg    13 6674.9 6741.5 -3324.5   6648.9                    
## intmodel3_neg   13 6679.7 6746.3 -3326.8   6653.7     0  0